{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np  # 矩阵操作\n",
    "import pandas as pd # SQL数据处理\n",
    "from sklearn.metrics import r2_score  #评价回归预测模型的性能\n",
    "import matplotlib.pyplot as plt   #画图\n",
    "import seaborn as sns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = pd.read_csv('day.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>instant</th>\n",
       "      <th>dteday</th>\n",
       "      <th>season</th>\n",
       "      <th>yr</th>\n",
       "      <th>mnth</th>\n",
       "      <th>holiday</th>\n",
       "      <th>weekday</th>\n",
       "      <th>workingday</th>\n",
       "      <th>weathersit</th>\n",
       "      <th>temp</th>\n",
       "      <th>atemp</th>\n",
       "      <th>hum</th>\n",
       "      <th>windspeed</th>\n",
       "      <th>casual</th>\n",
       "      <th>registered</th>\n",
       "      <th>cnt</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>2011-01-01</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>6</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>0.344167</td>\n",
       "      <td>0.363625</td>\n",
       "      <td>0.805833</td>\n",
       "      <td>0.160446</td>\n",
       "      <td>331</td>\n",
       "      <td>654</td>\n",
       "      <td>985</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>2011-01-02</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>0.363478</td>\n",
       "      <td>0.353739</td>\n",
       "      <td>0.696087</td>\n",
       "      <td>0.248539</td>\n",
       "      <td>131</td>\n",
       "      <td>670</td>\n",
       "      <td>801</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>2011-01-03</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0.196364</td>\n",
       "      <td>0.189405</td>\n",
       "      <td>0.437273</td>\n",
       "      <td>0.248309</td>\n",
       "      <td>120</td>\n",
       "      <td>1229</td>\n",
       "      <td>1349</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>4</td>\n",
       "      <td>2011-01-04</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0.200000</td>\n",
       "      <td>0.212122</td>\n",
       "      <td>0.590435</td>\n",
       "      <td>0.160296</td>\n",
       "      <td>108</td>\n",
       "      <td>1454</td>\n",
       "      <td>1562</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>5</td>\n",
       "      <td>2011-01-05</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>3</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0.226957</td>\n",
       "      <td>0.229270</td>\n",
       "      <td>0.436957</td>\n",
       "      <td>0.186900</td>\n",
       "      <td>82</td>\n",
       "      <td>1518</td>\n",
       "      <td>1600</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   instant      dteday  season  yr  mnth  holiday  weekday  workingday  \\\n",
       "0        1  2011-01-01       1   0     1        0        6           0   \n",
       "1        2  2011-01-02       1   0     1        0        0           0   \n",
       "2        3  2011-01-03       1   0     1        0        1           1   \n",
       "3        4  2011-01-04       1   0     1        0        2           1   \n",
       "4        5  2011-01-05       1   0     1        0        3           1   \n",
       "\n",
       "   weathersit      temp     atemp       hum  windspeed  casual  registered  \\\n",
       "0           2  0.344167  0.363625  0.805833   0.160446     331         654   \n",
       "1           2  0.363478  0.353739  0.696087   0.248539     131         670   \n",
       "2           1  0.196364  0.189405  0.437273   0.248309     120        1229   \n",
       "3           1  0.200000  0.212122  0.590435   0.160296     108        1454   \n",
       "4           1  0.226957  0.229270  0.436957   0.186900      82        1518   \n",
       "\n",
       "    cnt  \n",
       "0   985  \n",
       "1   801  \n",
       "2  1349  \n",
       "3  1562  \n",
       "4  1600  "
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(731, 16)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>instant</th>\n",
       "      <th>season</th>\n",
       "      <th>yr</th>\n",
       "      <th>mnth</th>\n",
       "      <th>holiday</th>\n",
       "      <th>weekday</th>\n",
       "      <th>workingday</th>\n",
       "      <th>weathersit</th>\n",
       "      <th>temp</th>\n",
       "      <th>atemp</th>\n",
       "      <th>hum</th>\n",
       "      <th>windspeed</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>6</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>0.344167</td>\n",
       "      <td>0.363625</td>\n",
       "      <td>0.805833</td>\n",
       "      <td>0.160446</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>0.363478</td>\n",
       "      <td>0.353739</td>\n",
       "      <td>0.696087</td>\n",
       "      <td>0.248539</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0.196364</td>\n",
       "      <td>0.189405</td>\n",
       "      <td>0.437273</td>\n",
       "      <td>0.248309</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>4</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0.200000</td>\n",
       "      <td>0.212122</td>\n",
       "      <td>0.590435</td>\n",
       "      <td>0.160296</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>5</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>3</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0.226957</td>\n",
       "      <td>0.229270</td>\n",
       "      <td>0.436957</td>\n",
       "      <td>0.186900</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   instant  season  yr  mnth  holiday  weekday  workingday  weathersit  \\\n",
       "0        1       1   0     1        0        6           0           2   \n",
       "1        2       1   0     1        0        0           0           2   \n",
       "2        3       1   0     1        0        1           1           1   \n",
       "3        4       1   0     1        0        2           1           1   \n",
       "4        5       1   0     1        0        3           1           1   \n",
       "\n",
       "       temp     atemp       hum  windspeed  \n",
       "0  0.344167  0.363625  0.805833   0.160446  \n",
       "1  0.363478  0.353739  0.696087   0.248539  \n",
       "2  0.196364  0.189405  0.437273   0.248309  \n",
       "3  0.200000  0.212122  0.590435   0.160296  \n",
       "4  0.226957  0.229270  0.436957   0.186900  "
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "y = data['cnt'].values\n",
    "X = data.drop('cnt', axis = 1)\n",
    "X = X.drop('dteday', axis = 1)\n",
    "X = X.drop('casual', axis = 1)\n",
    "X = X.drop('registered', axis = 1)   ###去掉原数据中的cnt，dteday，casual，registered这几列\n",
    "X.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(584, 12)\n",
      "(584,)\n",
      "(147, 12)\n",
      "(147,)\n"
     ]
    }
   ],
   "source": [
    "from sklearn.model_selection import train_test_split\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=33, test_size=0.2)\n",
    "print (X_train.shape)\n",
    "print (y_train.shape)\n",
    "print (X_test.shape)\n",
    "print (y_test.shape)    ##讲原数据分为训练和测试"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/anaconda2/lib/python2.7/site-packages/sklearn/utils/validation.py:475: DataConversionWarning: Data with input dtype int64 was converted to float64 by StandardScaler.\n",
      "  warnings.warn(msg, DataConversionWarning)\n",
      "/anaconda2/lib/python2.7/site-packages/sklearn/utils/validation.py:475: DataConversionWarning: Data with input dtype int64 was converted to float64 by StandardScaler.\n",
      "  warnings.warn(msg, DataConversionWarning)\n",
      "/anaconda2/lib/python2.7/site-packages/sklearn/utils/validation.py:475: DataConversionWarning: Data with input dtype int64 was converted to float64 by StandardScaler.\n",
      "  warnings.warn(msg, DataConversionWarning)\n"
     ]
    }
   ],
   "source": [
    "from sklearn.preprocessing import StandardScaler\n",
    "\n",
    "# 分别初始化对特征和目标值的标准化器\n",
    "ss_X = StandardScaler()\n",
    "ss_y = StandardScaler()\n",
    "\n",
    "# 分别对训练和测试数据的特征以及目标值进行标准化处理\n",
    "X_train = ss_X.fit_transform(X_train)\n",
    "X_test = ss_X.transform(X_test)\n",
    "\n",
    "#对y做标准化不是必须\n",
    "#对y标准化的好处是不同问题的w差异不太大，同时正则参数的范围也有限\n",
    "y_train = ss_y.fit_transform(y_train.reshape(-1, 1))\n",
    "y_test = ss_y.transform(y_test.reshape(-1, 1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[-2.94481124e-16]\n",
      "[[-1.04327382  0.30419853  1.4286816   0.42743455 -0.02917687  0.06956207\n",
      "   0.04514525 -0.19855368  0.23569164  0.24877962 -0.05353419 -0.0891602 ]]\n"
     ]
    }
   ],
   "source": [
    "from sklearn.linear_model import LinearRegression\n",
    "linreg = LinearRegression()\n",
    "linreg.fit(X_train, y_train)\n",
    "y_test_pred_lr = linreg.predict(X_test)\n",
    "y_train_pred_lr = linreg.predict(X_train)\n",
    "print (linreg.intercept_)\n",
    "print (linreg.coef_)          ##得到coefficients"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "('MSE:', 0.18268928443777904)\n",
      "('RMSE:', 0.4274216705289743)\n"
     ]
    }
   ],
   "source": [
    "y_pred = linreg.predict(X_test)\n",
    "from sklearn import metrics\n",
    "\n",
    "print (\"MSE:\",metrics.mean_squared_error(y_test, y_pred))\n",
    "\n",
    "print (\"RMSE:\",np.sqrt(metrics.mean_squared_error(y_test, y_pred)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfAAAAEYCAYAAACju6QJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJztnXl8VOXVx39nliSTkA2ILAMIaAuCLJEgIBaFoiggjaCixVbburXVKq2pWK1AbRXFpe2rb99irRvUgqARRMUFqIpFARNEEBBRwAQlCMOSDMlk5nn/mNzhzsxdZ+6sOd/Pxw9mcuc+596Z3POcnYQQYBiGYRgms7ClWgCGYRiGYczDCpxhGIZhMhBW4AzDMAyTgbACZxiGYZgMhBU4wzAMw2QgrMAZhmEYJgNhBc4wDMMwGQgrcIZhGIbJQFiBMwzDMEwG4ki1AGbo3Lmz6N27d6rFYBiGYZiEsWnTpoNCiDK94zJKgffu3RsbN25MtRgMwzAMkzCIaI+R49iFzjAMwzAZCCtwhmEYhslAWIEzDMMwTAbCCpxhGIZhMhBW4AzDMAyTgaRMgRNRTyJaQ0SfEtFWIro1VbIwDMMwTKaRyjKyVgC/EUJ8RESFADYR0ZtCiG0plIlhGIZhMoKUWeBCiP1CiI/a/v8YgE8BuFMlD8MwDMNkEmkRAyei3gDKAXyg8LsbiGgjEW1saGhItmgMwzBMFrJlyxY8/vjjqRYjLlKuwImoA4BlAG4TQhyN/L0QYoEQokIIUVFWpttZjmGYNKa6pg6j561Gn1krMXrealTX1KVaJKadsWfPHlxzzTUYMmQIZs+ejaNHo9ROxpBSBU5ETgSV9yIhxIuplIVhmMRSXVOHO1/cgjqPFwJAnceLO1/c0m6UOG9eUsuhQ4cwc+ZMfPe738XixYtx++23Y+fOnSgqKkq1aDGTyix0AvAkgE+FEI+kSg6GYZLD/FU74PX5w17z+vyYv2pHiiRKHu1985IONDY2YsGCBbj66qvx2Wef4cEHH0THjh1TLVZcpNICHw3gRwDGEVFt238TUygPwzAJpN7jNfV6NtGeNy+pwufz4X//938xY8YMAEDPnj2xd+9ePPnkk+jZs2eKpbOGlJWRCSHeA0CpWp9hmOTSvcSFOgVl3b3EldB1q2vqMH/VDtR7vOhe4kLVhH6oLE9uwYvaJqXO40V1TV3S5clmAoEAlixZgrvvvhuff/45xowZg2PHjqGwsBCdOnVKtXiWkvIkNoZh2gdVE/rB5bSHveZy2lE1oV/C1kwX17XWJoVd6dbx2WefoaKiAldddRUKCgqwcuVKrF27FoWFhakWLSGwAmcYJilUlrtx/9RBcJe4QADcJS7cP3VQQq3PdHFdK21eUilPtnH8+HEAQNeuXZGTk4PnnnsONTU1mDhxIoLpVtlJKjuxMQzTzqgsdyfVXZyouLtZt7z0u9sW1yZEnvbKzp07cdddd+GTTz7Bli1bUFhYiPXr16darKTBFjjDMFmLmus6nrh7rG75ynI33AmQpz1SX1+PG2+8EQMGDMDrr7+O6dOno7W1NdViJR1W4AzDZC2JiLvH45ZPRR5AtlFTU4PTTz8dTz31FH7xi1/g888/x5w5c5CXl5dq0ZIOu9AZhslaJNe1lVno8bjlEyFPe8Dr9WLr1q2oqKjA4MGDccstt+DGG29E3759Uy1aSiEhRKplMExFRYXYuHFjqsVgGCaL0Ytvj563WrEczl3iwrpZ45IpalqUyCWS1tZWPPPMM5gzZw6ampqwd+9eFBQUpFqshENEm4QQFXrHsQudYRimDSPx7XRxg6dLiVwiEELgpZdewuDBg3Hdddehe/fuWLZsWbtQ3mZgC5xhGEvJZKtQzboucTlRkOsIXdPY/mVYs70hpdcYiycgUz6bDRs24Oyzz0a/fv1w33334dJLL83qcrBIjFrgHANnGMYyJKtQSvKSrEIAaakoIlGLY3u8Pni8PgDBa1q2qS7hNex6mI3FW/HZJHID8PHHH+PDDz/Eddddh+HDh2PlypW48MIL4XCwmlKDXegMw1hGujROiRWj5VypvCZpqpma71TtGuL9bBLlsv/iiy/wox/9CEOHDsVdd92FpqYmAMDEiRNZeevACpxhGMvI9IElWh3TIonlmuIdKSpXokpExuLl66m9x+h1WL05O3jwIG699Vb069cPS5cuRVVVFbZv3478/PyYztce4e0NwzCWkaqBJVahVObV1NKKw02+qGPNXpMVLmwlJSrhjnBpR66nhtHrsHpz5vF4sGDBAlx77bW455570KNHj5jO055hBc4wjGVUTegXpTQyrVFJZLtXJUUYyzVpWbBGFbiasiQgKnFNS9lLmLmOeDdnLS0tWLBgATZv3ownnngCp59+Ovbt24fOnTsber9RMiVRzwrYhc4wjGWkYmBJorHqmqywYM20htU6byzXEWv5XCAQwKJFi9C/f3/ccsst2LVrF7zeoGyJUN7ZWlqnBFvgDMNYSrIHlliFluVmxTVZEV4w4+FQWy/WhjNK4YWx/cswf9UOzFxcq2jtbtu2DT/84Q+xefNmDBkyBK+99homTJiQsJIwK7wcmQQrcIZh2j3JKH+zIrxgphXr2P5lWLh+b9TrnqYWVNfUxXRd8o2M1j274LslKCgoQLdu3eBwOLBo0SLk9fse7lz5KW5a+yqAYG39nCkDLVWsmZ5EaRZW4AzDtHuSYblZ1QfdqDdgzfYGxdcbW/yWbE6U7tnRr7/EjCvvBR07gGdXrMbUYT2xYcMGvFxbj6qlm+Hznyx+83h9qHphc9xyyMn0JEqzsAJnGKbdkyzLLZnhBS3ZtTYnWqEE+e/kdeitxw7iyHv/wvEtb4GcuSg6eyruXFoDm82GynI35q/aEaa8JXwBYekmScnLQQh6I7IRVuAMw7R7stFyU7smCSUFr+UWB6BYltZctx3f/Pt3ECKAwrMmo/ic6bDnF6NZIKSctTYTVm6SKsvd2LjnEBat3xvaYAgAyzbVoeLUjlkXB+csdIZh2j3pMqDESvSa0pTkO6Ne0wolyH8X8DWj5ZvPAQA5XU9Hh/KJcF//d3QcfwPs+cWh90rKWWsjZPUmac32hqgudZnUDdAMrMAZhmn3ZFv5m+Tq9vr8UMv3Pn6iNaq8SiuUUO/xQgT8OFb7OuoXXI8DL8yBaG0B2R3oOO46OIq7RL1PUs5qLmwbYPkmqT0lsrELnWEYBplb/hZJpBtcrWe6UvxZze3erTgPh7a+h89f+wdaD32F3O79UXL+tSBHjqoccg+GWkJdntOmWYYWC9kYDlGDLXCGYZgswkgHNolIq7RqQj847eE2u9NOmFTmwY5Fc2Cz2VA29W50uXo+8nqeqXreSA+GmvXb5AuENV2pWro57qYr2RgOUYMtcIZhmCzCjKtYKQ4umewt3+xGS8MXKB08HgPOGoGXX34Zvu5D8chbu6Ky0CORlOXoeatR7/HCRgS/0HpHEJ9fYO6KrXFZ4VaV62UCrMAZhmEsIF16cOtln8uRdKoke53HC5/na3jefQ5N2/4De4dOKOg/Bg+9sRPrZk0BAEyr6AUAGDr3jdCM9EjmLN+K5tZAyBNgRHlLHG7yGWo0k+jOeZkAK3CGYZg4SUYnN6Mo1UKrccTrC8l+3PMtjrz/bxyrfR1ks6No5OUoHjEN5HCizuMNWdOSstTqhqqm2O1ECAihu8mQl64pKel0ut+phISJnVGqqaioEBs3bky1GAzDZBDJsIxHz1utqJBK853Iz3FYsraZ61BruBJJvtOG5lYBvxBoafgS+5+5DR0GjUfxOVfBUdhJ9X0up91wnD2SL+dNAqBtwUsQwpPwXE477p86KOQtiCTWPu/pBhFtEkJU6B7HCpxhmGxFbRSo1SVifWat1FSU8a4dz3WobS5Eqw/Hal9F6+H96HjBTQAAf6MH9oISQzLZVeLaBIAICKjckD9PHxqyoqte2Ayf2oEquEtcqhsTAvBF2wYhkzGqwDkLnWGYrEWrMYmVGC1RinXteK4jMitbBPw4/slq1P3jJhx++wn4vv0Kwh+0hI0qbyAY1470ohOAGSN7qSpvAKHxnpXlbsy/fAjsJieTSR4IJbKxVEwLVuAMw2QtyWrqodf1LN6147kOqUkNgKCb/Olb8e3KR2DP64BTrrgXp0y/F2QPz0Y3ci2R7m1JeVec2lG1eQwQvvGoLHfj4SuGGL53AELhg/ZSKqYFK3CGYbKWZFlqSp3cSlwKJVoxrq32HhuRobrpC/uVwl3iCrY5JRs6T/ktul7zKFx9yqNmc9uJQteiRqTyRtvPa7Y3YP6qHbrhBPnGQ37v9JCUdLZ1zosVzkJnGCZrsWIGt1EiS5fU4taxrK2WWe4XQjP7etu2bfjd736HAwcOoOqxJfjdSy2wXfuXKKUtQQCuGtEzdC6lNV1OG7y+gOL7jXoXtDYxSpsDILixmDYsKJc8I/7Rtph6e4QVOMMwWUsqm3pYubb0nt8s2RyVOKY0GnTfvn2YPXs2nnnmGRQUFOC3v/0tJg/qAiIKk6d3Jxfe//yQ5uQu+fFj+5dh2SZ1i19SzHp16I3NJ/uwz12xFYebTmajq1nvfiGw+MN9WLxhX2g0qVr5mFWVB+lS268GZ6EzDMNkCGrZ7lL2dXVNHe7+2xJs/ecdIACTr7wW//zzfejcubPi+dQy1NXKsdSOB4LehWnD3Fj58f4whayG00YAQXFOuFnk8lpVeZCsCgYljGahswXOMAyTIag1QOmST/jLkrfwf1ta0VTUG4VDL0bR8Ep81rkb3tvXjEpl/W06OU7LRW4jYPGH+wyXhZktH9NC3mhGqW2rkpdCD63M/3SxwlmBMwzDpIBY3LORsXDhb0Xz1rexc8MS3NEKdLnu7yBHDjqOvwFAUOH8Zslm1WlfZid3aXVQa2yJrbGLFRBOuu3V2raazf7X28Skg3uds9AZhmGSjOSelU/ikuqjtZCyr7sX56Fp+3s48NQv8c1r/4MzvnsaSif9GmSPtsn8QqiuYbYcy0y5XCJw2ihqWppa0lskapuS6po6jJ63Gn1mrcToeatD90drExPr52c1rMAZhmFUUHu4x3tsLI1ZpPPPXFwLz2cb0fDyPJzWtQQvv/wy3nvvPfQ9UzdkGrWG2XIsMyVfVuC0E0pczpBs8y8fgvmXDQmT12gHPKVNiZYi1trcJKtBkB4pTWIjon8CmAzggBBCfbhsG5zExjBMsjCTxGTkWCP9ydVagVbX1GHm4y/h+IG9KBhwHoQQ8O/+AH+54/rQdDAlGcysYRathDY91KzmfKcNpQW5ptzSanLIB6eonUcviU/NTa6XTBgvmZLE9jSAxwA8m2I5GIZhwjCTxKR3rFHlquS23bVrF6679np8+/Fa2Iu7IL//uSCbHY7TRuKRt3aFFLgk05zlWzWHhFjVxMbM1LNI1DYwXl8A20wOI1Gr9TeSLa4X51YbS2o2dyBRpNSFLoR4B8ChVMrAMAwTSXVNnap1qfTQ11MESgo+kkg37zfffINf/vKXOOOMM3Do0/+iaNR0dP/JX0G2k25dpXWbW5WbrAAnk730XPxayF35NnNtzHWJRQHG05Ut1k596dLKNe1j4ER0AxFtJKKNDQ0NqRaHYZgsR7KW1VB6uOspAr0M6NJ8Z5TS2bt3L5544glcf/31GHb7cygd8yPYcgvC3ieAMGWstVGQu62Vkq6MxPAjY8ZWZp477aSqAPVkqyx3Y92scXh0+lAAwMzFtYY2KbEq4nRp5ZpqF7ouQogFABYAwRh4isVhGCbL0VKCag93vZatWuVXAJCf48BFZ3TCo48+ivr6esyfPx/Dhw/Hvn370KVLl2BS1dLNik1P5N3ItDYKke/UcvGrdTgz4kmIlYIch6ICNCqb0ePkxNMtT829nkzSXoEzTKaTDvWijHG0lKCalaWnCLTixSLgx453X0HxA1ehxXMAF110EZZt2INH3v48rIWpVrq1VO9t1sLRcvErxfutnuIm54hK3N6obGrHzV2xVfPvLx0UcaywAmeYBBKLVcBYQ6wbJy1rWT4GMxItRaDWy7z56134duWj8B3cg5yup6PnD2/DsMsn4+7ln4Z9Zxat36urnNUamLicduQ5bYrtTfVc/JGv63kS4kEtDGFUNrXjDjf5QteebX9/KY2BE9HzAP4LoB8RfUVEP0ulPAxjNelSL9reiKfRhlazklgbdkibCb8QIAABXzMAtI33JHSecge6/vgR2HoOxvMf7Iv6zsQTO8xz2jBpcDfVWG91TR1sKtPJIpVqohq5aMWdjSaaGU2Ay6a/v1RnoV8lhOgmhHAKIXoIIZ5MpTwMYzWx9JRm4ieejZNesxKvz4/bDCZJAeGbiZaGL/HNsj+gYdlcCCHgKCpDt5/8DwrO+B6Igo9jNUs6Vg43+bB4wz7ItwFS0hwQHBmqtKaSUlVK3tJCS8FIWwa9BDCjiWZmNhfZ8vfHLnSGSSDpUi/a3oh34yS5w9UadgDG3bHzV+3AsW/3w/PuIjR+shqU40LxyMtgh0AAFDWb264wjEMLI61EfX4RlgB3om2et1pSmp1IM94vf12roUsAwfnhJ3wBFLucIAI8TT7TyWKSrFrhEKXjGptbFWvis+XvjxU4w1hIZNxVmp+slp2ciDU5Sc66jZNezNfIdKrPa97HN8v+AAAoGl6JolGXw+4qalNu9qjvxrRhbsXvjI2Uy7byc+wICJjKDpfkVlW8Qqhek5HvuJwTvgAenT40ru+k0USzyOPUOuQlu147UaR9HTjDZApKcddlm+owbZg7YfWi6TJUId2Ip9GGvOa4qaU1OLdaAyWrvrGxEdu2bQMA9B5QjsIhE+C+4e8oHfcz2F1FAE5+FyK/G3+sHKT4epNKzXVjiz/sO2a0uUqdxwu1Q6WBHZG111rfcTUEEBW6MNM3Ph7SpV47UaS0F7pZuBc6k87o9VXOljUzhVg8E0oWm9NOKMhxqLYnld/rpR9+iTvu+zP2vPUsnK4CLHxtHex2e8ytPuVouarlMgyd+4ZmK1UJNVc9AZgxspeiFyDXYVM8txQLV5NP3iPcTI/59kqm9EJnmKwhmQlrknIy0+4zk4lFGcdS36sUE/b5BQpyHZg8pFtUOZdk1QcCAfz2wb/jrw/+Eb7D9cjtMRDF512LWS9+gvmXD8H9Uwcpym/0uqpr6tDY3Koqd53Hiz6zVqJ7icuQ8o503csRAFZ+vF8xCVDtPfUeLx6dPhQzF9cqxuPloQszPeYZbViBM4xFJCthzchgjGxJ0gGSW0uvtvGRXMVy5UQApg0LbhKqq6vx8J2/gLPzqSi7bDZcfStARPAFBG5bXAt3W6x4zfYG1Hu8mL9qBzbuORRm5RrtMKaG5NLWS2ojBBvSzF2xVbE2vDTfqfi6Ft1LXKgsd2PjnkOqmxwJrsywDo6BM4xFJGvAgV47y2xK0gGSW0uvtvGxE0XJcOLrXVi69EUAwCWXXILOlXei20/+ivzThkdlltd5vFi4fm9Y7HjR+r2K1/WbJZvDYsNm25cKQDW2DQAOO2HjnkM4fiLaonfaCVpR1dJ8p+Z3/I+Vg/Do9KGhmHOJy4k8py2sN3msA0SYaNgCZxiLiKevshm0LBV3ErPQk5X9brXFpiW3Wk9z+c++Q3XwvPMcmna8h4MdeyAQ+D3sdjsK+o02JYeanpTi0pJFHuu4TreKR8jnF3j+g32K8e+CHIdqS1MAmH3JQADa33EpdKHmOVHLsk/2pjMbqjdYgTOMhSSjr7Kaqz6ZiWuJdmvLH642lWSrWCw2LbmBk9a+lOAlbYjmr9qBPV/V4ci6f+H45jdAjhwUn3Ml+l3wQ9hsQUdmLK5nPeSyRKLlKpe+C2p17Gp15ke8PhS7nIpxdJfTFvpsjXzGap6TNdsbVHMCkkW2tDhmBc4wGYbe5KtkkMhEpMiHq9EuYUbQGnhxwhcIW1NaQ7qeWx/dhOMfv4nC8otRfM6VcBSU4oKhvULnmX3JQNWJYfEgtV+NjL+rrWJmClok3UtcaGpRTpbLM9lCVctzYmSjm0gLOVsS6ViBM0yGkSxXvRaJTETSi/kW5Njxp0sHRSV6GbkfWgMvImnyenHbXfdi07Ay3Hvvvdg49SI8U9QV9g4dAQQV6LJNdag4tWOYQpKqA4x0SJNDBNX4sxTX1nKNS8jLscb2L8PC9XsNrS8p/pmLaxV/74khsS3WpM5EW8jZkkjHCpxhMpBUj0BMZMa93kO0sSXYi3z+qh0Y278Mr2zeH+by1XrYG7FIRcCPxk9Ww/Pev+A/1oDNjksQCASwZntDSHlLRFpt8s9FvqkwoshLXM4wL0CUXDjpGteq/5df85rtDZpr2okQECJs06NWnmj2s43HU5RoCzlbWhxzFjrDMGGodcnS61DmtBMam1vj7q5l9CEqZXYrxWvVstTH9i/TzNBurt+B/f+8Bd++9hfYO5Ri4HUPYfny5bDZbKattspyN6om9DN8PZ4mn+YQFflaRise9DZDASHwxbxJWDdrXFginxXVFPF0QUu0hZysipFEwxY4w7Rz5JZiscuJxpbWUBy3zuPFzMW1uG1xbZhL+HCTD3Ybhb3m84uQMo3H5alkucVCpIVVXVMXVcstIVpbQI4c2FyFwfGelXei08Dv4Y/TBoeOMWu1Ga3flp9HsuDVLGxpLaNhlBKdxDol2a0M0cTqKUq0hZwOYSgrYAXOMO2YSCWjZM2KiH8l/AFtx7DX58ec5VujHoryLnKR2d7yB/6c5VsNdRVTg9rWkj+sI5VpS8OX8PznGYAIp0y7B87S7uj208dARJg3bXCY7EZdwtU1dbqyR8bH5edR67oWuZYR5ahV061lcaY6RJOMRM1UX6MVsAJnmHaM2SYhZvF4fSFXulJyV2TNM3DywVpZ7sbd1VsMJ2FFIg3RkB7Scvdr65Fv4Hl3IRq3roUtNx9FIy+HEAJEJ8d7Si54+UM+12EL3S8bhbvqpdrnqhc2w6exuSEAj04fqtpaVclqL813YvYlA00rHK2a7nTuPZ4tFnKi0VXgRNQFwH0AugshLiaiAQBGCSGeTLh0DNOOSEVjiWRk3c5ZvhXNrSeTs9RUm9d3MjlNuvaKUzti8YZ9MZdm1Xm8IStccss27foQDdX3gciGohFTUTTyctjzOii+V14jHqlYJR0tP27+qh2ayhsId5VHMmf5VsUNVX6Ow1JXdGSyWzqSDRZyojFigT8N4CkAd7X9vBPAYgCswBnGIlLVWMJsnXAsmHWDRylEFeWt1uAkkjtf3AJvUyMGFzai3mNDbo8BKBwyAUUjLoejqLPme+UWtpanQtp86OFy2jG2fxlGz1utaH2r3atYN1rp0DOASRxGstA7CyGWAAgAgBCiFUDifG4M0w5JZr9vOWYf5FIGd2m+0/Dc6ViQmqtobS6uGtFTM6McAITfhwMfvIwfTxiJp+69DQEhYM/rgI4X/FxXeUvUe7wWeioEFm/Ypzi/XeuzjjV5K9vnYbd3jFjgjUTUCW2eLyIaCeBIQqVimHZGqhpLVJa7VadSRRLZZ93o3OlYOdzkg41Ouqrl2AiqGeUAIEQATZ++A8+7C9Hq+Rq5Pc9E6XnXRg0ZMYKkPK3wVHh9AYXX/KHQiRrxWMzsis5ejCjwXwNYDuA0IloHoAzAZQmVimHaGbGUKFkVL599yUDdcielPutaCVJWoRZODghtl3bTp+/i4IqH4CzrjVMum4O8vsNiUt6Ek8rTitI2Neo8XtUOa6X5TlbAjCK6ClwI8RERnQegH4Lf5x1CiMT/5TJMO8JMrNJsvFxP2eu1AI2UQzqftR2/46e5fgf8TR7knz4C+f1Go7PNjvx+54DIpmrJa0EAZozsFXav5PXyPn8AjS3KCl1vmIgSjc2tcNopLObvctpDE8AYJhISOkkgRPRLAIuEEJ62n0sBXCWE+N8kyBdGRUWF2LhxY7KXZZikYNSq1mqjGWklK5UluZx2zThoZGMXomCXsO4lLoztXxY1CjKR5DttEKAo+XMdtpD73vftV/C88yyadr4PZ1lvdPvJ/8RkbcvRK9vSatIiv79qn5UaThuhQ54jdL8TXYmQDSM1sxEi2iSEqNA9zoACrxVCDI14rUYIUR6njKZhBc5kE7E+PNWsOgLwxbxJYa+ZUfZK8kUqKbMDOpQwmj0OnBzwIW/4IvU/P3jg6+B4z4/fBDlzUTT8UhQNr4QtN1/3fEZkjOwRLkftvtqJ8PAVQ8L6oZt1vSdrLGwsmzsmORhV4EZi4DYiItGm6YnIDiAnXgEZJh1IlQUST9mYmXh5PMlxSpnx8SrvP08P2gJGlZqkbKXxnr07ubBo/V4IAL4DX+D4lrdReNYkFI+aDntBie4Gw+C+QbXBjITa/QsIoRqeUGpVq0SyJmJly0jN9owRBb4KwBIi+j8E/zZuAvB6QqVimCSQqtprIPaHZ3VNneK8ZrV4eaw9patr6iyvDy9xORXjyd1LXGhsbtXNaG9sasKr/1oCBPwoHnUF8voOg/umf8BReLIcrCTfiUmDu+H5D/YZtvL1UPpczNzXyCxweStZJWxEYS1gE0W2jNRszxipA78DwGoAPwfwSwBvA/htIoVimGSQqtprILaHp7ThiCz5KnE5Vd2esUxdktaxEpfTjjlT1JOxJg/pFiWnhAj4cWzzKtQvuAGetU+h5etdobancuUNBEvPFq7fa5nyloj8XOKZZiVNKXPaleP0fiFCteGJRG0Tl2kjNdszRrLQAwD+1vYfw2QNqbRAYrGM1fqWF+Q6FC08ybqdNsyNNdsbQsNDIvt3G10nHuQbDCXPx7JNdTirVzHe//xQmAv8xFef4tvX/oLWQ18hp1s/dL7kduT1GmR6fTsRCvMcMdetR34u8fbq1uowByTHlc1d2jIfVQVOREuEEFcQ0RYohJWEEIMV3sYwGUOiRxZqEcvD08iGQ005ThvmDsse1woXWL2Biey7reb5WL/78MnJZ60+kMMJW64LZLOh7NLfwfWdUTFnlweEwJwpAw21O41E7XOJp0GKkXucjCY+AA8MyWS0LPBb2/6dnAxBGEaOVnKZVYlnqbRAYnl4GtlwqClHpZiwmpWnN0PaDEr3U00x+YVAy4HdOLz2GdhcL3kCAAAgAElEQVRy81H2gzuQU9Yb3X76eNxlYdIAkZlLajUT2aQs9zXbGxKq1Iz0oE/GRpK7tGU2qgpcCLG/LeP8SSHE+CTKxLRztJLLAFiWeJYsC0Rtw2H24Vk1oR+qlm4Oc7067RSmILWUoxJKx8cbPpYywUvaashnRkwYU9og+Dxf4+h7C3F863+C4z1HXRE13jNWpE1EdU2d5rURkJTyLUB58yiHXdmMETRj4EIIPxE1EVGxEIL7nzNJQS+5zMrSl0RbIJZnukcqoIifzU4XU7LyYmmRKiltqV86EL3Rmrm4VtGF3bhjHQ4unw+y2VA0chqKRlymON7TKC6nDXlOe1gzFEkeLZKZvKVUXiZvmMOubMYIRsrITgDYQkRvAmiUXhRC/CphUjHtmliSy9K19MXKWlulWdO+gAg7l5Jlp1YbLe/zLSeWEaOPTh8adj2j563WrCEPNDfB33QEztJuyOsxEB2GTEDxqMujsspjwesLwOsLoMTlDClCJXnkRHoykgG7r5l4MaLAV7b9xzBJQS/Wm6rEs1iwMtPdyLmUwgJqylhA2QtQNaEfql7YHLVZUCMySQ1Qn9wlWn04Vvsajvx3MRwlXdH16odgLyhBpwt/bmgtM3i8vpDVrXW/9dqmMqmB27zqY6SM7BkiygHQH8G/+R1CiJaES8a0W/SSy9Kx9EXtYWM0093Iw8rouSItO61+3KPnrVYcbmJ0xGik5VpdU4e5K7ZGHSdEAI3b/oMj7y5E65FvkNtrMErPuybu+LYekrdD7d4lq20pY45UNlnKJHQbuRDRRACfA/grgMcA7CKiixMtGNN+qSx34/6pg+AucYEQfMhKdcRav0sV0sOmzuOFwMmHTXVNnaGGH1rvlzO2f5ni+mqvSyjJIKG2lsdgFvr04T2j6ruVFH/jlrfx7SsPg3ILcMrlc9Hlyj8ht7v+posQjGnHQ73HG1fjFSb5pLLJUiZhxIX+CICxQohdAEBEpyHoUn8tkYIx7Rut+GC6xQ61HjaSdadlXRuNk6/Z3qC4vtrrEpHjQiMx0ypUa+3I62iu245AcyNcfYehYMB5oBxXaLynEUpcTsyZMhAzY6jdluNy2kKyyYeisEs2feE2r8YwosAPSMq7jd0ADiRIHobJOPQeNnobDqMPq3geapIMvWcpp7PIlbVav3W190l9u6Vz+A7uw+F3n4V353+R0+27yOtzFsiRg4L+50a9X21Otw3AnCkDDc0d1xtg0uQLoKlNNqmcztPEUcB0JpVNljIJIwp8KxG9CmAJgn8nlwPYQERTAUAI8WIC5WOYtEftYWN0KIXRh5UVDzWtUZ6j562Oad73nS9uwcY9h+A/ehCedf/C8S1vgZy5KD53BoqGV6rGuUs1GsaItvpxPeXtjiFjHgAaW/yoWroZAMdU0xFu82oMI76sPADfADgPwPkAGgB0BHAJ4uzSRkQXEdEOItpFRLPiORfDpAq1GLPRoRRG47PxxHGra+owet5qzSEfdR4vFq3fq6q8nTZlRSx1emvevxPHt65G4bBL4L7xHygZfRVsOcqbCwKQn6NuPwihbVUTgqNJ180aB3eMVpnPLzimmqakY65LOkLC4qk9hhcOdnnbCeACAF8B2ADgKiHENrX3VFRUiI0bNyZJQoYxTnVNHX6zZLOigjSS6Wy0ZMZMaY18bKWem9kIkRZzwHcCxzatANkcKDr7Uggh4D/2LRxF+rXcV4/sFZrrHQtXj+yFP1YGh5rcXb0FC9fvVTxO77oJwBfzJsUoBcMkBiLaJISo0DvOiAs9UZwNYJcQYjcAENG/AfwAgKoCZ9oXmVQHWlnuVk22MhOjNnucZFlH3qPIMhwrtumHm3xwl7jw1aHjOP7xmziy7l/wHz+E/DPGAEBwvKcB5e1y2lBxasfQhDSzOG2EilM7hn5WS+KzE+GqET01NwocU2UymVQqcDeAfbKfvwIwIvIgIroBwA0A0KtXr+RIxqScTKwDTXbijdY9SsRIUAJwmn8PNj75e7QcqkNu9/7oPOW3yOt5pqnzeH0B3PnilqgJaUaJ7D6ntkEKCBGy0pUs9FR0XwMya2PKpDfxFVjGh1JATWls6QIhRIUQoqKsTLvelckeMrEOVC9GLVnLfWatxOh5q3Vj43po3SOry22E3wcBoLbeiy6lHdD/R39Al6vnm1becjnXbG8IxTkB5QeCGvLrU9sgSa//sXIQ/jx9KEpcztDvSvOdmH/ZkKQrTqM1/wxjBK154L/WeqMQ4pE41/4KQE/Zzz0A1Md5TiZLyMQ6UK3pZonwKGjdo1j6mSvFi1u++RyH1z4Ne2EndJ54G44Vnoovd26DzWbT7PBmhHqPNxQSMHsuudI2krEcT+8AKy1mK3vjM4yWC72w7d9+AIYDWN728yUA3rFg7Q0AvkNEfQDUAbgSwA8tOC+TBWiVZvWZtTJtXY9qiiLeB7c8IU2rFAwI3iOlxDW9hC7573yH98Pz7kI0ffof2PIKUdx3GIDgrPDvPbgW9R4vSvKdcNrIcM/0SLqXuMKuSw2X066rnIHEjIW1euOViRtTJn3Rmgc+FwCI6A0AZwkhjrX9PAfAC/EuLIRoJaKbAawCYAfwTyFEdBNlpl2iNi9ZUlyZEBOXE8+DO1KJaClv+e/lR7lLXBjbvwzPf7BP9/2N2/6DgysfCWaXj7oCxWdPhS2vA5x2wvETraFM9MNNPjjthBKXEx6vz1Smu8tpx9j+ZZozsSW5qyb001XOierOZ7XFzA1KGCsxksTWC4C8bVELgN5WLC6EeBXAq1aci8kuIq0qm4LVabXrMRZXqdH3xPPgjjchTRqrWbVUucyNAPibmxA4cRyO4lOQ22MgOgy5CMXnTIejw8ls74IcBzwRs8J9foGCXAcKch2GXeBypax1XZKlncrWuVZbzNyghLESIwr8OQAfEtFLCG6wLwXwbEKlYhiEW1V9VFqAWuV6jMVVauY98Ty4471Gj9eHuSu2wuePVt6i1YejNa/iyH8XI+eUPuhy5Z/gKOocNd7TXeKyRJlJiWR6XdbSpVe51RZzIt39TPvDyDjRPxHRawC+1/bST4QQNYkVi2HCSbTrMRZXqZn3qD24ASjWccuJJSEtksiWpSLgR+O2tfC8uwj+oweQd+pQlJx3jer7JYvZ7Jz2SDxeX5QVH0k6jfhMhMWcbsN4mMzFaB14PoCjQoiniKiMiPoIIb5IpGAMIyfRrsdYrEuz71FqwlK1dHPIMq7zeBX7c1dN6IeqFzbHnCymxPHa13Hozb8hp8tp6HTRLXD1Kdc8XpIn8jNw2ghNLcG4uBXd3mwENDa3oveslWkxOUxp4zW2fxnmr9qBmYtr2YJmUoquAiei2QAqEMxGfwqAE8BCAKMTK1pi4WYKmYXVrsfIz7+4LRErEi0LPx6vQHVNHWYuqUVkSNrnF5i7YmvYdVWWuzF3xVbVwR9GOfHVp4Dfh7xTB6PgzO/Dll9saLynVKcd+RkUu5xobDmZ1KbXstSIcg8IhD6HdElYlG+8MrHBEJO96PZCJ6JaAOUAPhJClLe99rEQYnAS5AvDql7okX+EQNCa42b57QOlz99pJ0AgzMqVfyeUNnxAtEWq9T0y05ucgLBNSp9ZK1WPJwCkMpYTAFoa9sDzzrPw7voAuT3PRNcfzlNdV6lkS+164q0DN0s6uNbVrjkdZGOyB6O90I10YmsRQS0v2k5cEK9wqSYTu3wx1qH0+fv8Ah3yHIrTj9S6ZwHA/VMHoTT/ZIevXIfyn5T8HIC+NRrZpUvNqneXuPDFvEkoynNG/a71aAMOvvpn7H/qFpzYuwUl3/sRTrlsjuqapflOUxOgkl27nA610lzHzaQTRmLgS4jo7wBKiOh6AD8F8I/EipVY+I+wfaP2OXuafKi558Kw19SmjEkbvqoJ/XDCFzh5Dq8PMxfX4rbFtSjNd0II4IjXp1gGZwT5Olo5AEru/xN7t6Bx21oUVkxB8agrYHcVqa7jtBNmXzLQVIJVicY870SQDrXSXMfNpBNGstAfIqILABxFMA5+jxDizYRLlkD4jzB7iCWXwejnL1nNaoq33uNVtOalo+XKLRblLV9HLweAAPhbTuDYpuWw5XVAYflEFAw4D3m9BsFRpDxDoCDHjsYWf9v/hz8KjNzXWC5Jr4ucGgSkRa0013Ez6YSRJLYHhBB3AHhT4bWMhP8Is4NYE4qMfv56jUa6a9RGxwKRslLsLksiU7qupR9+GazlXvc8/I2HUXDmeBSWTwTZ7KrK22kLj5l7vL7QvQMQdV9vW1yLuSu2hqx0IOhZ0LwehIcKXE675gQyAnD6KQXYdaAxqgXsjJG90iI/heu4mXTCiAv9AgCRyvpihdcyBv4jzA5ibXNp9PPXUs6Swtfr422UP08fCkA5KU5rY7l69Wr86Iof48S3dcjtMQCdK+9EXo8Buuv5AoAvoHzvGptbFRXs4SZf2AZJzZMhJXQpWfFqmyI7ER6+YohqwmA6/W1yHTeTLmhNI/s5gF8AOI2IPpb9qhDA+4kWLNHwH2HmE08ug5HPX32gSjBZbebiWhS7nHDaSbHLmVHcJa6QLBv3HAr1K7cTYdowZTl9Ph+cTiccDgf85EDZtHvgOm04iMwM5YxGbzMi3yCp9auvP+JF71krw+q3767eophLIBEQInSd/LfJMMbQykL/F4KTx15u+1f6b5gQYkYSZGMYTfTmQEdidh630nxvp51gJ4LHG5yP7fH64u5eIp8XvmxTXUjJ+YXAsk11YXJ+9NFHuPDCCzFz5kwAwJgxY1Ax8x/IP/3sKOVtI2D0aR2hRI5dWdHbDOh/aYNUWe6OysIHToYBpJDGjCf+i4Xr92rGvjn/hGHMo6rAhRBHhBBfAvgLgENCiD1CiD0AfEQ0IlkCMowaSgpWzeWsVgqmpcQlBSUvqyrIcUR1RPMFgtZyLLictjCXvlJI4DdLNqPHjU+g85CxGDZsGD766CP0798/dMy4AV0Uz/3DEb3w5bfKFrVaVzcjzd7kyray3I38HPVInNfnx7rPD2mej/NPGCY2jMTA/wbgLNnPjQqvMe2YVMUszeQyxBMvl/9ebaiKX4ioJihGkJegqbn+j3z8Jr59/X9Adgc6nXsV/vKn32PGmDNCv1/58X7F963Z3qB6zjiS4qOUbTyJfOkytIRhMhEjCpyErF2bECJAREZ7qDNZTqpbSxqNl6opmTqPV3eYiBytxC15UptULlWa78TxE62qFq/cmpWfO9DciECzF46izsjtMRCFQy9C8agrYe9Qivlr9oUUeHVNnWottnRNSvKqlXOVuJw40hYeUKLE5bRs2IqdiLuXMUwcGOnEtpuIfkVEzrb/bgWwO9GCMcYwG9e1Giu62iXjGtRirASYcqtrue0ry91YN2scvpw3CZ/fPxFfzpuEmnsuxPSze0LJwx7pOq6a0A+55MfRD19C3f9dh0NvPA4AcJZ2Q8cLfg57h1IAwbj70LlvhDwfWtesJu9VI3oqvj5nykDNkP6cKQMN3RMjXDWip+n3MAxzEiMK/CYA5wCoA/AVgBEAbkikUIwxYonrWk28Xe2SdQ1qSiZSWeltPpTi4tOGuTF/1Q7FDYiUmBZp7Ja4nGFtSv1+Pzyb38Thp3+Bw2ueRE7X01F8rnquqFS3rWX5SpsKpfaof6wcpNo21a2y2bFRcI535DXK1wBgKB/g6pG98MfKQbrHMQyjju4wk3TCqmEm2UI6DFaIV4ZkXoPRYSIE4It5kwyfU2ugidHre/jhh3H77bejoqICdadPRc6pxmYFabnCa2dfqPAOfZSuSYkSlxNzpgxUDTnoXbuVuRPpXjvOMGaIe5gJEf227d//IaK/Rv5npbBMbKRDT3czmeBKJPMaJBe3u8Sl6Sa2ERl25+uFELSu77333sO6desAAD/72c/wwgsv4MMPPzSsvIGTyXNyJFd4rERa7WoWtcfrw22La1H+hzcU75PWd8NKz0s6eKIYJhVoJaN92vYvm7xpSjr0dI+3q52Za7DKytLbHBidQ11dU6fqwpbWUJoz3tLwJU68vxDfe2A9Sr4zHCVTZ7ddzyi8XFtv6lrkyXNWWp/y5EC1zHuJyA5t8nMAyt+N0fNWx1QVoESsFQYMk+moKnAhxIq2f59JnjiMGazs6R6Pcoync5bRa7Ay291M1rSaIpDk0VqjuqYOjS2toddajxyA572FaPxkDXLzO6Bs7E+QVz4pZDXOXFxrqieMPHkukYrKyP1Su0+RStyIZ8Is6eCJYphUoNVKdQU0ekwJIaYkRCLGMFb1dI9VOVphERu9hnitLLmsSu1PtWq4lRSB1qATeZ90+RreLz5C46fvouycaej9/Rk40BLRwUznGv48fajp+23FZ6TWMjUSpfuk9t1S8kwAsXmP0sETxTCpQMuF/lDbv1MBdAWwsO3nqwB8mUCZGBNYYX3FohyttIjl1yApnJmLa8MUTjxWVqSsHq8PThuhNN8JT5MvbNCGUUWgte79U4PZ1fsOHMLRDdVwFHZGh8EXoMOg8XD1HQZnURkaWnTFDkPql27m3lr1GUnHzl2xVXP+t9J9Uvtu5TltUZumWL1HPF2Qaa9oudD/AwBEdK8QYozsVyuI6J2ES8YkjViUYyLijkoKp+qFzZi7YquqdWrEylKS1RcQyM9xoOae8Exto4pAzeorcTnR6vPh5tkPoeHdRQg0etBhyEXoMPgCkN0BR1FZSGajbvxIGYxa1VZ+RtLmobqmDnOWb42yntXuk9p3yNPkw6MxeBTUZAN4uiDT/jDSUa2MiPoKIXYDABH1AaA8ZJjJSGJxQSYi7qimaLWsvsbmVlTX1Gk+rLW6sPWZtTLqgW9EEVRN6IeqFzZHdVhr2P4hZvzlWrQc3o/cnmei9NK7kOs+2fZUruiMuKVdThvynMHJZ/NX7cDY/mVh87S1rGozn5HRTYH0mlyJl+Y7w+aEy9H6blkZu+cJZkx7xIgCnwlgLRFJ3dd6A7gxYRIxSScWF2RJvlNRscYTd4xF+UsNTQB1t7BWEpa87Eg6hxFFUFnuDrmUhRBAwA+yO9DqD0A4cnHKZbOR17ciakKYvHkLAN26dK8vAG9bv/Q6jxeL1u9VbT4TKbfaZ1QSMT3MjKtdqUZc3s89EnZvM0zi0FXgQojXieg7AKTxR9uFEM2JFYtJJmZdkNU1dTh+ojXqdaed4nowx9pTW88tbCQJS167LY/1ajUr8TT50Lx/Jzz/eRo5XU5H6difIq/PWejWpxxEyi0WJEtafp8IQaXqaVLvQS6h9nulzY9aj6bI18242s265dm9zTCJQ1eBE1E+gF8DOFUIcT0RfYeI+gkhXkm8eEyyMOOCnL9qh+JwjoIcR1wPZqPZzkpoWe9Kbl8l6jxeVC3dHJY57vH6UPXC5rDzAMDOnTtxbOWD+PaTd2BzFSG/37kA0GZxq7cSlSz+qhc2A4TQWlphAiMoeT6OqFxr5OtmXO2xhE7Yvc0wicFIL/SnALQAGNX281cA/pgwiZi0RzUxyeuLayBJZE9tMxhx3Te3qrt6gaDalStvCV9AhPVH//vf/44BAwag8fON6PS9H8J94z9QWD4RQNA9XOJyRp1D6ZxKa8WCNJAl8r6r3RMbkaHjlF43cyzDMInFiAI/TQjxIAAfAAghvNAyMZisR+thHW8rS3m7UyWcNuUvX1NLq+Z6WnXbwfOSpvt639cN2L9/P6pr6vB/253IH3Ixzpz5DG769Sx0Ki0OHZfrsGHykG5w2pL3JyLJHXnf1Qa4+IXQPU4tTh1v61yGYazDiAJvISIX2p4TRHQaAI6Bt2OMjI80O1I0EjUr3xdQjgMfbvKhaulmDJ37RpQXQKvlKRCsse6QpxxNCviaceSDF7F/wfW4/JobceeLW3A4tws6XnATGvwuLFy/F0dl+QAerw8L1+9Fq8r8b7NcPbKXqePl913yaCj1Mlc6TmkyWSRmjmUYJrEYyUKfDeB1AD2JaBGA0QCuTaRQTHoTmZhkJrHKKLEktPn8IhTjlqzRjXsOYdkmdctcmowV2e9bBPxo/ORteN77F/zHDqL8nPNxvP9kRSver6CsrVDfdiJUnNoRr2zerxm7j0R+3yvL3Zi5uNbQccloncswjHVoWuAUzMjZjmA3tmsBPA+gQgixNuGSMYpU19Rh9LzVccWarUBydX8xb5KquzueuKgRK18Pr8+P5z/Yp+o6JwBj+wdbGkTKemT9C/j2tb/C3qEj7l2wBB+tW4Mj+dYqrdJ8Z+jeKTncJVd3S6u5pL7Ia1H7HIoNxOoZhklfNBW4CA4LrxZCfCuEWCmEeEUIcTBJsjERpOvYxETEReNJaJOjNCtbQgBYtqkO1TV1qJrQD6J+G5r37wQAFA69GO7L7sa/X3kbd19/OQBrE7VcTjtmXzIQ62aNw5fzJuHR6UNVXd1NGnXWkRAQdd+rJvRTjMl7vD70TvFGkGGY2DESA19PRMMTLgmji97saT0SZb0nKi4qWfnxpIOpzbKW8Pr8mPP0q3ji7huw97nfonnDUhCAXt274rHf3YRLz+oROtaqRC07UdT9qSx3I6Cx2VAi8soIwIyRvRQngqnF+IHYN4Lp4g1imPaKkRj4WAA3EdGXABoRfE4IIcTgRArGRLe31Js9rXcuq4aPKJHIuKjatdsIKMpz4ojXh5J8J46faA2rT3c57Zg2zB3WelRO65Fv4Hl3IfZsXYs9JcWYN28ebrnlFuTn5yvKUVnuxm0q8WQlnDYKq/WWZFLb3KhdZ2m+Eyd8gahuZtOGubFme4OhBikenTpzsz3SE/19YhhGHyMK/OKES8FEofSAVGu3qefara6pw2+WbI5yJ8c7fEQ6t9HOZbGi1OCFAAQEUJDrCK2n1s+74tSOipPGmnb+F0071qH7mCvwSfXfUFpaqiuLW0XJumUTzeTrA8a7kKm1HZ19yUBT51HCSFKgmaTDRAyzYRjGHFrzwPMA3ATgdABbADwphIjun8kkBKUHpACilLherFnaCKjFguPJFK+uqTPcuSwe5FnvkhKS1z5XLT25nlrpU2W5G8+v24lf/e5eoLgbCgach8Lyieh45hg8dO04Q8ob0O7trbW+2etUUtSJ7nJnJsafiGE2DMOYQysG/gyACgSV98UAHrZqUSK6nIi2ElGAiCqsOm82ofYgFICpWLNeA5N4ErPmr9phqHOZFUjx8NL86Mxpn19g7oqtqu9taWnBY489htumfg8H31kIx8FdIAA9OhfjoWvHmZ6Nncg6aHl2/7pZ5mTTO688KTAyfm426ZA7sjFM6tFyoQ8QQgwCACJ6EsCHFq77CYKlaX+38JxZhZrLU6pbNoqWRRRvprjWuWOxxIyMtFTrGX64yYfR81ZHvefVV1/FLbfcgt27d+P888/H8uXLMWLECNOyyTET7zc6pjMZyOWOVy6eMsYwqUdLgYeelEKI1sixiPEghPgUQNSoReYkVj0g1TYCSpnQZh/qWnFVs5aYFUlR0nuEELhkcFc4HA40NzejsLAQr732GiZMmJDU71w6J3rFm3TIU8YYJvWQUImNEpEfwaxzIOhxcwFowsks9KK4FydaC+B2IcRGI8dXVFSIjRsNHZoVWGG9Kc1vVsqENnpc5LkjY+BAMPt6/uVDdGWVX5+NSDFOH+lxGDr3Dc2uZM31O9C07ln8+seXYvbs2RBCQAgBm81IxaS1jJ632hIvCsMw7Qsi2iSE0A0vq1rgQoi42mAR0VsAuir86i4hxMsmznMDgBsAoFcvc32hMx0rSrOMWkqxZBVLr8eShR65YTCaZDdnykBUvbA5apyp79t98LzzHJp2vg9bfgnc7uD6RJQyT0+qE73SyX3PMIz1GCkjiwkhxHiLzrMAwAIgaIFbcc72hpGNgJpSqfN40WfWSlUFEOsmQy+5TiLSFa+UkX500wocfvsJkDMXxefOQP/xV+K66yablslq1EIMyUj0Smf3PcMw1pB8vyKTliRyRKgSRqxQtZh/Zbkbr9xYjjkX9IDLaUeu+wwUnjUJ7hueQLfzr8asKeWWyBgvqRy9GW/XPoZh0p+UKHAiupSIvgIwCsBKIlqVCjkykUS1rzQ6InTOcvVyLTPkObW/emrlWV6vFw888AD69u2Ltc89ivunDkLf/oPQafyN6OXuphuzT2brz1SO3ky1+55hmMSTMBe6FkKIlwC8lIq1M5lEukWNjgj1eH2orqmLq4wKALwaAzrkjVEkWltb8dRTT2HOnDmor6/HxIkT8etf/xqDBxtz4afKpZyq0ZupdN8zDJMcVLPQ05H2loUeSTKzmtXWMrOeWmZ7ntOmWs+ttsZdd92F++67DyNHjsQDDzyAMWPG6K4t3zgcamxW3DRI62RbwlcsVQUMw6QHcWehM+lHMt2iVRP6qQ7uMLqeWhzWSPJavceLtWvXorS0FEOGDMEvfvELDB8+HD/4wQ9ARJoKV8na1lonGxO+uE6bYbIfVuApIFZrL5lu0cpyd1h5mJn1pOvTG56hRss3u+F55xmMfWATOg0+H/94eiEqy92h0jA9hWs0w126lmwdzJEq9z3DMMmBs9CTjKR86trizGayu5Od1Tz7koGK643tX6aaDCa/PrP4PF/j4IqHsP/pX6G5fgdKzv8p8sffEnV/9DKszXgkqib044QvhmEyErbAk0w81l6i3KJqHgGl9cb2Lwubry1tQDbuOYQ12xtitroBoGnbf9C0878oGnk5ikdMgy2vA4Dw+1NdU6c7F93I6EwgOGdbstg54YthmEyDFXiSidfas9otqueOjhyAoTZXfNH6vaqZ62oEmptwdMNLyDmlL/K/OwqFFT9AwaDxcBR2ijpWHqtWQ1K4RkZnyuds82AOhmEyEVbgSSbdynv0PALyeHbkLHI5ZpS3aPXhWO1rOPLfxQg0HUHR2VOR/91RsOXkwZaTp/getVi1BAEY278MgLKnYmz/MqzZ3qA5Z5sTvjiBUhcAABB3SURBVBiGySS4jCzJxFreY3WZk16iGQF4dPpQXUvWLE2ffYDDby9A65FvkHfqYJScdy1yu31X8z1OO2H+ZUMwc3Gt5kaBy6QYhskGuIwsTYnF2rO6zElpExEFQbWMzAwEICAEIAIgmx0O/wl06lSKTlN/jeNlA+Cw2VQHmUgU5Dg0Y9US2ZA5zjAMYxTOQk8BleVurJs1Dl/Mm4R1s8bpKhyr+1obKbMy45ix26KnfTlthD9PH4qbzvDD88JdOLahGu4SFx77/S2o/2wr/nTLDOTnOHSVNwAcaRsfaqTdK2eOMwzTXmALPAOwuszJaiXnD0QrYcexejzzhydQXV2NLl26YN41Y3DttSc7q5mt1QaUJ5GpHcswDJPtsALPAKxOfDNaZhUrRz54EXv+8zTqOhTg3nvvxW233YYOHTqEfq9VChaJ005h2eBSVrxaLgFnjjMM015gF3oGoOQ6lmddW3G+WHA57ShxOQEAfu8x+L1HAQC57v7ods6l2L17N+6+++4o5a1VCkYyb3xpvhPzLxuiGGJI5aQvhmGYdICz0DOEu6u3RNVaO22EDnkOeJp8pjPTjZaHhdayE6YP7xlVitXs9eKWu/+Eg++/gIIB56PThT/XzAbXGpIiyeHmMi6GYdoxnIWeZsRbBrZme0OUkvUFRKhXuZHMdDUZ9ErKnDZEWcI+nw///Oc/MXfuXDTs34/SM0bBVX6xrvLVir9L15cNw0QYhmESDbvQk0A8/c8ljCSeeX1+zF2xVbFPuZYMUla8WyWmfkqRK0qR3nHHHbjpppvQp08fvPvuuzi07X08fnMlAGDm4tqoHukSRuP2Wln21TV1qr3YGYZh2gtsgScBK6ZdGU08O9zkU7TKjcigl+2+evVqdO3aFQMGDMDNN9+M888/H5dccklovGfV0s3w+UVo7aqlm0M90tV6qWuhJE82jv5kGIaJBbbAk4AVZWCxJp5JStqIDMVtCWmROA9/iQkTJuD73/8+HnzwQQBA3759MWXKFFBb1tncFVtDylvC5xdYuH5vmNW/bFMdpg1zhyWfleYrr6tkrVtdE88wDJOpsAWeBKwoA4vs4AYY7z8uWb96MlBEPxbf4Xp43nkOTdvfxbcdO+Lhhx9Gj3N+gNHzVkfF0ZXmhivh9fnx/Af7EBAi9H4AhkvCePQnwzBMELbALUAvJmvVHG95BzczSIpSqRStzuPFaXe+it6zVkYp4eNb3oL38w9RPGo6du/ejb5jp2P2ys/iiuUDgF+IsPcDMFwSprbp4QYuDMO0N9gCjxMjMVm9/ufVNXWYu2JrSIGWuJyYM2WgZkxXzaJWKglramnFxj2HkOuwhVm50nFSO9NAcxOOfrAMuT0GwNV3GIpHTEPhWZNxag83iouLMX/VJlX3dYnLCY/XmBWu9H6llrJKWfNmR39aPQSGYRgmXeA68DhRq2t2l7iwbtY4hXeEE5n8JeG0EeZfrtzERHqfkiKbNsyNVzbvN6VMRasPx2peDY739B5F0agrUDrmxyE5pp8drP/Wm1xW9cJm+BTaqupBQJRXQWtqG2BsGEysk98YhmFSCdeBJ4l4Y7LzV+2IUt5AsMZbK0s90qovyXdCCGDR+r2wRQazNWja8T4Orf4H/EcPIO/UoSg57xrkdvtO6PetAYHFH+7TVMzdS1xh8pht02o2Wc3IABi9c7ACZxgm02EFHifxJqhpKXq9TYBaX3C9CV9Br4sAkQ1+71HY84vQ6eJfwdV7aPSxgKbylruvJXn6zFqpmmAX6eJPZLIaJ7wxDJPNcBJbnMSboKal6I1uAsxM9jrx1af45l934HjNqwCADoMvQNcfP6KovPVwl7gwbVhwTrc8gU9NbgIwY2SvpCWrccIbwzDZDFvgcaKXoKZH1YR+qjFwo5sAIxZlS8MeeN55Ft5dH8BeUAqnqwOcNsCH2IaaSC1TlRL4pg1zRzVrkZT3HysHGTq/2WS1RJ2DYRgmXWEFbgGS6zjW9wIwnYUuR82NbydCQAgEal7E1289DZvThZIxP0a/71+B8YN7aSamSdgoeB65G11Sgmox5jXbG3D/1EGheLidCH4hsGZ7Q6h1qx7xboysOgfDMEy6wlnoGYTWMJJISzPH14jZPzgTM8YMwFtvvYXXX38dd955Jzp16qSand2jNA+fHWgMvZbrsOGBaYMBKCtBtVi3lFXOWeAMwzDm4Sz0LMNIvfn8VTvw1YHDwNZXUff+C6i13YgZY+Zj/PjxGD9+fOhcapbzLpnyBhDKZlfzMOgl8HEWOMMwTOLgJLYMQa8H+KQzT8HVxTvQ8vzN2PvGPzF+3Fhce+21iudSi5lHWtN6Pcb1Evg4C5xhGCZxsAWeJOLtCKanDG+99Vb87W9/w7nnnoulS5di9OjRqK6pww1tfculOvEjXh9sbTHpeNYF9GPMVvSAZxiGYZRhBZ4ElNzfty2uxZzlWw0nqykpQ++XtXD36AUA+NWvfoWJEydi0qRJofGe8jXlfc6VlLdSC1ZpXS20Evg4C5xhGCZxsAs9CajVaXu8PsPDQOTu6uavd+Gbf9+NA4vvxil73gQA9O/fH5MnTw6N9zRSG24nCtVjzxjZy5KBK3Iqy92Gh5QwDMMw5mALPAlouaGNJnVVlrtRv3c37vn9Pfh2y1o48ovw09/Mwf/+aZbpNSUCQoT1IK84taPlJVfxlNgxDMMw6rACTwJqsWAJo0ldn619ESd2b8Dvf/973H777SgqKop5TekYOaxsGYZhMgdW4ElAKRYsRy3OfOTIETz44INw9SnHawc7Yl9gBM64dTTOunSUpvI2sqYR9ziP4mQYhklfWIEnAaVuaxJKivTEiRN4/PHHcd999+HQoUPoPOZqFIy6EjZXIRr8iKr/1lozclrZEa/PkDI2UnfOMAzDpA7uxJZk9KzaJUuW4Pbbb8e+fftw4YUX4uvvXoojBT2jzmN03nisxDvnnGEYhokN7sSWpijFmYUQEELAZrOhvr4eXbt2xdNPP41x48ahz6yViudJdDMUbsLCMAyT3qSkjIyI5hPRdiL6mIheIqKSVMiRDrz77rs499xz8eyzzwIAbr75ZnzwwQcYNy5o5SZrJGZ1TR1Gz1sdGgta7HImZV2GYRgmNlJVB/4mgDOFEIMB7ARwZ4rkSBlbtmzB5MmTMWbMGHzxxRdwuYKK0eFwhGq5gfjnjRtBinfXebwQCMa7G1ta4bRR2HHchIVhGCZ9SIkCF0K8IYRobftxPYAeqZAjVcyZMwdDhgzBe++9h/vvvx+7du3C9OnTFY9NRjMUpaYvPr9AhzwHN2FhGIZJU9IhBv5TAIvVfklENwC4AQB69eqVLJkSSkVFBW6//XbMmjULHTt21D0+0fXZanFtT5MPNfdcmLB1GYZhmNhJmAInorcAdFX41V1CiJfbjrkLQCuARWrnEUIsALAACGahJ0DUpDN58mRMnjzZknNZUavNQ0cYhmEyj4QpcCHEeK3fE9E1ACYD+L7IpFq2NMKqWm0eOsIwDJN5pCoL/SIAdwCYIoRoSoUM2YDejHCj8NARhmGYzCNVMfDHAOQCeLMt43q9EOKmFMmSsVhZq8190BmGYTKLlChwIcTpqVg32+DYNcMwTPuF54FnMMmoEWcYhmHSE1bgGURktzQAHLtmGIZpp6RDHThjALWM8/unDuLhIgzDMO0QtsAzBKsyzhmGYZjsgBV4hsDTwRiGYRg5rMAzhGRNJWMYhmEyA1bgGQJnnDMMwzByOIktQ5Ayy+Pte84wDMNkB6zAMwjulsYwDMNIsAudYRiGYTIQVuAMwzAMk4GwC91irJjPzTAMwzB6sAK3EKvmczMMwzCMHuxCtxDulsYwDMMkC1bgFsLd0hiGYZhkwQrcQrhbGsMwDJMsWIFbCHdLYxiGYZIFJ7FZCHdLYxiGYZIFK3CL4W5pDMMwTDJgFzrDMAzDZCCswBmGYRgmA2EFzjAMwzAZCCtwhmEYhslAWIEzDMMwTAbCCpxhGIZhMhASQqRaBsMQUQOAPamWwyI6AziYaiGyHL7HiYfvceLhe5x40u0enyqEKNM7KKMUeDZBRBuFEBWpliOb4XucePgeJx6+x4knU+8xu9AZhmEYJgNhBc4wDMMwGQgr8NSxINUCtAP4HicevseJh+9x4snIe8wxcIZhGIbJQNgCZxiGYZgMhBU4wzAMw2QgrMBTBBHNJ6LtRPQxEb1ERCWplinbIKLLiWgrEQWIKONKRNIZIrqIiHYQ0S4impVqebIRIvonER0gok9SLUs2QkQ9iWgNEX3a9py4NdUymYUVeOp4E8CZQojBAHYCuDPF8mQjnwCYCuCdVAuSTRCRHcDjAC4GMADAVUQ0ILVSZSVPA7go1UJkMa0AfiOEOAPASAC/zLTvMSvwFCGEeEMI0dr243oAPVIpTzYihPhUCLEj1XJkIWcD2CWE2C2EaAHwbwA/SLFMWYcQ4h0Ah1ItR7YihNgvhPio7f+PAfgUgDu1UpmDFXh68FMAr6VaCIYxiBvAPtnPXyHDHnwMI4eIegMoB/BBaiUxhyPVAmQzRPQWgK4Kv7pLCPFy2zF3IejKWZRM2bIFI/eYsRxSeI3rUZmMhIg6AFgG4DYhxNFUy2MGVuAJRAgxXuv3RHQNgMkAvi+4ID8m9O4xkxC+AtBT9nMPAPUpkoVhYoaInAgq70VCiBdTLY9Z2IWeIojoIgB3AJgihGhKtTwMY4INAL5DRH2IKAfAlQCWp1gmhjEFERGAJwF8KoR4JNXyxAIr8NTxGIBCAG8SUS0R/V+qBco2iOhSIvoKwCgAK4loVaplygbaki9vBrAKwcSfJUKIramVKvsgoucB/BdAPyL6ioh+lmqZsozRAH4EYFzbM7iWiCamWigzcCtVhmEYhslA2AJnGIZhmAyEFTjDMAzDZCCswBmGYRgmA2EFzjAMwzAZCCtwhmEYhslAuJELw2QoRNQJwNttP3YF4AfQ0Pbz2W19ypMt0yoAl7X1lmYYJoFwGRnDZAFENAfAcSHEQxGvE4J/54EEr5+UdRiGOQm70BkmyyCi04nok7bmQB8B6ElEHtnvrySif7T9fxciepGINhLRh0Q0UuF817XNrF/VNgP8bpV1urU1HClp+/1P2ubdbyaip4yuxzCMMdiFzjDZyQAAPxFC3EREWn/nfwXwoBBifdtEplcAnKlw3Nltr7cA2EBErwA4Ll8HAIKGOEBEQxBsFXyOEOIQEXU0uR7DMDqwAmeY7ORzIcQGA8eNR7BVp/RzKRG5hBDeiONWCSEOAwARVQM4F8DrGuuMA7BYCHEIAKR/TazHMIwOrMAZJjtplP1/AOEjQPNk/08wlvAWmSwj/dwYeaDsvEoJNkbXYxhGB46BM0yW05ZYdpiIvkNENgCXyn79FoBfSj8Q0VCV01xIRCVElA/gBwDW6Sz7FoArJde5zIVudD2GYXRgBc4w7YM7EHR5v43gPG+JXwIY3ZZstg3A9Srvfw/AvwDUAHheCFGrtZgQ4mMADwJ4h4hqAcw3uR7DMDpwGRnDMJoQ0XUAzhRC3JZqWRiGOQlb4AzDMAyTgbAFzjAMwzAZCFvgDMMwDJOBsAJnGIZhmAyEFTjDMAzDZCCswBmGYRgmA2EFzjAMwzAZyP8Dh8Xj5W0qAV4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 504x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(7, 4))\n",
    "plt.scatter(y_train, y_train_pred_lr)\n",
    "plt.plot([-2, 2], [-2, 2], '--k')   #数据已经标准化，3倍标准差即可\n",
    "plt.axis('tight')\n",
    "plt.xlabel('True price')\n",
    "plt.ylabel('Predicted price')\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "#岭回归／L2正则\n",
    "#class sklearn.linear_model.RidgeCV(alphas=(0.1, 1.0, 10.0), fit_intercept=True, \n",
    "#                                  normalize=False, scoring=None, cv=None, gcv_mode=None, \n",
    "#                                  store_cv_values=False)\n",
    "from sklearn.linear_model import  RidgeCV\n",
    "\n",
    "#设置超参数（正则参数）范围\n",
    "alphas = [ 0.01, 0.1, 1, 10,100]\n",
    "#n_alphas = 20\n",
    "#alphas = np.logspace(-5,2,n_alphas)\n",
    "\n",
    "#生成一个RidgeCV实例\n",
    "ridge = RidgeCV(alphas=alphas, store_cv_values=True)  \n",
    "\n",
    "#模型训练\n",
    "ridge.fit(X_train, y_train)    \n",
    "\n",
    "#预测\n",
    "y_test_pred_ridge = ridge.predict(X_test)\n",
    "y_train_pred_ridge = ridge.predict(X_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "('MSE:', 0.18104231340692375)\n",
      "('RMSE:', 0.4254906737014617)\n"
     ]
    }
   ],
   "source": [
    "\n",
    "print (\"MSE:\",metrics.mean_squared_error(y_test, y_test_pred_ridge))\n",
    "print (\"RMSE:\",np.sqrt(metrics.mean_squared_error(y_test, y_test_pred_ridge)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEKCAYAAAA4t9PUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3Xl8VfWd//HXh4SwhDUQtoRVNkH2iNqqiIJ7XVlsa6utjt38Tau/mam/sdPO2N9sddZ2utkZp9PaGUlALSquiLW1orlhk0VWxXsJEDaBBCEk+cwf94TepjfhxuRuyfv5ePDw3HO+95zPPULeOdvnmrsjIiLyUXVJdwEiIpLdFCQiItImChIREWkTBYmIiLSJgkRERNpEQSIiIm2iIBERkTZRkIiISJsoSEREpE1y011AKgwcONBHjRqV7jJERLJKRUXFQXcvPNu4ThEko0aNIhQKpbsMEZGsYma7ExmnU1siItImChIREWkTBYmIiLSJgkRERNpEQSIiIm2iIBERkTZRkIiISJsoSEREOqDIkRP8/fPvUHX8ZNK3pSAREemAykIRfvSrnZyu96RvS0EiItLBNDQ4SysiXDx2IEX9eiR9ewoSEZEO5vWdB9nzwYcsKhmeku0pSEREOpjSUIR+Pbty5eTBKdleUoPEzK42s61mtsPMHoiz/H4z22xmG8xspZmNDOZPN7M3zGxTsGxxzHvMzP7azLaZ2RYz++NkfgYRkWzywYlaXti0j5umF9EtNycl20xa918zywG+D8wHIkC5mS13980xw9YCJe5+wsy+BHwHWAycAD7r7tvNbBhQYWYvuPsHwJ3AcGCiuzeY2aBkfQYRkWzzy3WV1NY1pOy0FiT3iGQ2sMPdd7l7LfA4cGPsAHdf5e4ngpergeJg/jZ33x5MVwJVQGNP/C8BD7l7Q7C8KomfQUQkqywpD3NeUR8mDeuTsm0mM0iKgHDM60gwrzl3Ac81nWlms4E8YGcw6xxgsZmFzOw5MxvXTvWKiGS1jXuOsnnvMRan8GgEkhskFmde3Buazex2oAR4uMn8ocDPgc81HoEA3YCT7l4C/AR4tJl13hOETejAgQMf8SOIiGSP0lCYvNwu3DCtpd/Z218ygyRC9FpGo2KgsukgM5sHPAjc4O6nYub3AZ4FvuHuq5usd1kw/SQwNd7G3f0Rdy9x95LCwrN+U6SISFY7ebqep9bu4ZrzhtC3Z9eUbjuZQVIOjDOz0WaWB9wGLI8dYGYzgB8TDZGqmPl5REPiZ+5e1mS9TwGXB9NzgG1Jql9EJGu8sGkfx07WpfQie6Ok3bXl7nVmdi/wApADPOrum8zsISDk7suJnsrqBZSZGcD77n4DsAi4FBhgZncGq7zT3dcBfwf8wszuA6qBu5P1GUREskVpKExx/x5cNGZAyredtCABcPcVwIom874ZMz2vmfc9BjzWzLIPgOvasUwRkawWPnyC13cc4r554+nSJd7l6eTSk+0iIlmurCKCGSwoKU7L9hUkIiJZrL7BWRoKc8m4wpQ0aIxHQSIiksVe33GQyqMnWZSmoxFQkIiIZLXSUJh+Pbsyf1JqGjTGoyAREclSR2pqeXHT/pQ2aIxHQSIikqV+uW4PtfWpbdAYj4JERCQLuTtLQhGmFPVNaYPGeBQkIiJZaOOeY2zZe4xF56f3aAQUJCIiWak0FKZbbhdumDYs3aUoSEREss3J0/U8tS5o0NgjtQ0a41GQiIhkmRc27eN4mho0xqMgERHJMkvKwwwv6MGFaWjQGI+CREQki4QPn+C3Ow+xcNbwtDRojEdBIiKSRcpC4WiDxlnpa4nSlIJERCRL1Dc4SysiXDqukGFpatAYj4JERCRL/OZMg8bMuMjeSEEiIpIlSkNh+vfsyrxJg9Jdyu9RkIiIZIEjNbW8tGk/N81Ib4PGeBQkIiJZ4KmgQePiDGiJ0pSCREQkw7k7S8rDTC3uy8Qh6W3QGI+CREQkw7295yjv7DuecRfZGylIREQyXGODxk9kQIPGeBQkIiIZ7OTpen65rpJrpwzNiAaN8ShIREQy2PMbow0aF5ZkzpPsTSlIREQy2JLyMCMKenLh6Mxo0BiPgkREJEO9f+gEb+w6xKKS4oxp0BiPgkREJEOVVYTpYnBrBjVojEdBIiKSgc40aBxfyNC+mdOgMR4FiYhIBvr19gPszcAGjfEoSEREMlBZKEJBfh7zzh2c7lLOSkEiIpJhDtfU8uLmfdw0vYi83Mz/MZ35FYqIdDJPrt3D6XrPyAaN8ShIREQyiLtTFgozrbgvE4b0Tnc5CVGQiIhkkA2RoEFjlhyNgIJERCSjlIbCdO+auQ0a41GQiIhkiA9r61m+rpJrzxtKn+6Z2aAxHgWJiEiGeH7TXo6fqmNhFjw7EktBIiKSIZaUhxk5oCcXjilIdymtoiAREckAuw/VsHrXYRaVDMcscxs0xpPUIDGzq81sq5ntMLMH4iy/38w2m9kGM1tpZiOD+dPN7A0z2xQsWxznvd8zs+pk1i8ikiploUi0QePMzG7QGE/SgsTMcoDvA9cAk4BPmtmkJsPWAiXuPhVYCnwnmH8C+Ky7TwauBv7FzPrFrLsE6IeISAfQ2KBxzvhChvTtnu5yWi2ZRySzgR3uvsvda4HHgRtjB7j7Knc/EbxcDRQH87e5+/ZguhKoAgrhTEA9DPxZEmsXEUmZ17YfYN+x7GjQGE8yg6QICMe8jgTzmnMX8FzTmWY2G8gDdgaz7gWWu/veljZuZveYWcjMQgcOHGhV4SIiqVRaHqYgP48rsqBBYzzJDJJ4V4s87kCz24ESokcasfOHAj8HPufuDWY2DFgIfO9sG3f3R9y9xN1LCgsLW128iEgqHKo+xctb9nPzjOxo0BhPbhLXHQFij9OKgcqmg8xsHvAgMMfdT8XM7wM8C3zD3VcHs2cAY4EdwV0NPc1sh7uPTc5HEBFJrsYGjdl6WguSGyTlwDgzGw3sAW4DPhU7wMxmAD8Grnb3qpj5ecCTwM/cvaxxvrs/CwyJGVetEBGRbOXulIbCTBveL2saNMaTtOMod68jej3jBWALUOrum8zsITO7IRj2MNALKDOzdWa2PJi/CLgUuDOYv87MpierVhGRdFgfOcq2/dUszuKjEUjuEQnuvgJY0WTeN2Om5zXzvseAxxJYf6+21igiki6NDRqvnzY03aW0SXZe2RERyXIf1tbz9LpKrp2SXQ0a41GQiIikwXMbow0as/kieyMFiYhIGiwpDzNqQE8uGJ1dDRrjUZCIiKTYewdrePPdwyzMwgaN8ShIRERSrKwinLUNGuNRkIiIpFBjg8bLJgzKygaN8ShIRERS6LVtB9h/7BSLSjrG0QgoSEREUmpJeZgB+XlcPjE7GzTGoyAREUmRjtCgMZ6O80lERDLck2v3UNfgLDo/+58diaUgERFJAXdnSXmY6cP7MX5w9jZojEdBIiKSAuvCH7C9qprFHexoBBQkIiIpURqK0KNrDtdPze4GjfEoSEREkuxEbR1Pr482aOyd5Q0a41GQiIgk2XNv76P6VF2HPK0FChIRkaRbEgozemA+54/qn+5SkkJBIiKSRO8erOGtdw+zsKS4QzRojEdBIiKSRGWhjtWgMR4FiYhIktTVN7C0IsLcCYMY3KdjNGiMR0EiIpIkr20/QNXxUyzsAN+C2BIFiYhIkiwpDzOwVx5XnDso3aUklYJERCQJDlafYuWWKm6eUUTXnI79o7ZjfzoRkTR5ck3QoLGDn9YCBYmISLtzd0pDYWaM6Me4DtagMR4FiYhIO1vb2KCxExyNgIJERKTdlYXC9Oiaw3UdsEFjPAoSEZF2FG3QuJfrpnbMBo3xKEhERNrRig7eoDEeBYmISDsqLQ8zZmA+JSM7ZoPGeBQkIiLtZNeBat567zALS4Z32AaN8SQcJGZ2sZl9LpguNLPRyStLRCT7lFVEyOli3DqzKN2lpFRCQWJm3wK+Dvy/YFZX4LFkFSUikm3q6htYVhFh7oRCBnXgBo3xJHpEcjNwA1AD4O6VQMd/ykZEJEG/2tY5GjTGk2iQ1Lq7Aw5gZvnJK0lEJPs0Nmi8fGLHbtAYT6JBUmpmPwb6mdkfAS8DP0leWSIi2ePA8VO88k4Vt8ws7vANGuPJTWSQu/+Dmc0HjgETgG+6+0tJrUxEJEs8uTYSNGjsuN+C2JKEgiQ4lfWKu79kZhOACWbW1d1PJ7c8EZHMFm3QGGHmiH6MHdQ5Lx0negz2GtDNzIqIntb6HPDTZBUlIpIt1rz/ATuqqjvVk+xNJRok5u4ngFuA77n7zcCks77J7Goz22pmO8zsgTjL7zezzWa2wcxWmtnIYP50M3vDzDYFyxbHvOcXwTo3mtmjZtY5mtmISEYqC4XpmZfDdVOHpbuUtEk4SMzsIuDTwLPBvBZPi5lZDvB94BqiofNJM2saPmuBEnefCiwFvhPMPwF81t0nA1cD/2Jm/YJlvwAmAlOAHsDdCX4GEZF2VXOqjqfXV3LdlKH06pbQlYIOKdEg+SrwAPCEu28Knmp/5SzvmQ3scPdd7l4LPA7cGDvA3VcFRzoAq4HiYP42d98eTFcCVUBh8HqFB4C3Gt8jIpJqz769l5ra+k59WgsSvNhO9AihgehRxe2AETxT0oIiIBzzOgJc0ML4u4Dnms40s9lAHrCzyfyuwGeIhpyISMqVhcKMKcxnVidq0BhPokHyC+BPgI1EAyUR8TqWxQ2fIJxKgDlN5g8Ffg7c4e5Nt/sD4DV3/3Uz67wHuAdgxIgRCZYsIpKYnQeqKX/vCA9cM7FTNWiMJ9EgOeDuT7dy3REg9nivGKhsOsjM5gEPAnPc/VTM/D5Er8d8w91XN3nPt4ie6vpCcxt390eARwBKSkrOdvQkItIqZaFog8ZbOlmDxngSDZJvmdm/AyuBMz/s3f2JFt5TDowLrqfsAW4DPhU7wMxmAD8Grnb3qpj5ecCTwM/cvazJe+4GrgKuiHOUIiKSdHX1DSxbE2HuhEEM6t25GjTGk2iQfI7onVJd+d2pLQeaDRJ3rzOze4EXgBzg0eBC/UNAyN2XAw8DvYCy4NDwfXe/AVgEXAoMMLM7g1Xe6e7rgB8Bu4E3gvc84e4PJfg5RETa7NWtBzhw/FSnfZK9qUSDZJq7T2ntyt19BbCiybxvxkzPa+Z9j9FMm3p377z32IlIRlgSCjOwVzfmdsIGjfEkevvv6jjPgIiIdDpVx0/yyjtV3DqzqFM2aIwn0d/uLwbuMLN3iV4jMcCDBwlFRDqNJ9fsob7BO+X3jjQn0SC5OqlViIhkgWiDxjCzRvZn7KBe6S4nYyTaRn53sgsREcl0a94/ws4DNXzn1nPSXUpG0Qk+EZEElZZHggaNQ9NdSkZRkIiIJKDmVB3PbKjk+qlDye/EDRrjUZCIiCTg2Q1q0NgcBYmISAJKgwaNM0d07gaN8ShIRETOYkdVNaHdR1hcMrzTN2iMR0EiInIWZRVhcroYN6tBY1wKEhGRFpyub2BZxR4un6gGjc1RkIiItODVrQc4WH2KRXqSvVkKEhGRFiwpD1PYuxtzJxSmu5SMpSAREWlG1fGTrNpaxS0zi8hVg8Zmac+IiDTjiaBBo05rtUxBIiISR2ODxpKR/TmnUA0aW6IgERGJo2L3EXYdqGGRnmQ/KwWJiEgcS8rD5OflcN0UNWg8GwWJiEgT1afqePbtvVw/dZgaNCZAQSIi0sSzGyo5UVuv01oJUpCIiDRRGopwTmE+M0f0S3cpWUFBIiISY0fVcSp2H2Hx+WrQmCgFiYhIjLJQhNwuxs0zitNdStZQkIiIBE7XN7BsTYTLJw6isHe3dJeTNRQkIiKBVe9UcbC6Vk+yt5KCREQkUBqKNmi8TA0aW0VBIiICVB07yaqtB7h1ZrEaNLaS9paICLDsTINGXWRvLQWJiHR67k5ZKMzsUQWMUYPGVlOQiEinF9p9hF0Ha1ioo5GPRE1kWvDEmgh7j55kQH4eBfl5DOiVx4D8bhT0yqN3t1w9rCTSQZxp0DhVDRo/CgVJC57ZsJdX3qmKuywvpwsFvxcweRTkd4uZVvCIZIPqU3U8u2EvN04fRs88/Uj8KLTXWvDonedz8nQ9h2pqOVxdy8GaUxyuruVwTdPpWt47VMPh6lpqauvjrqul4CnIj76OTkfnKXhEUuOZ9ZV8eFoNGttCQXIW3bvmUNSvB0X9eiQ0PpHgOZRg8PTP78qA3wsbBY9IeysNhRk7qBczhqtB40elIGlnbQmeQzWnOBQTNoeqTyUcPF1zLE7Y/G66ID+PgUHwFOTn0ae7gkdkR9Vx1rz/AQ9ee67+PbSBgiTNPkrwHK6p5VAQPL+bruVwEESHamrZfegEh6pPnTV4CvK7BQHTGDbdFDzSaZQ2NmicWZTuUrKagiTLdO+aw7B+PRiWhOA5XFNL9am6uOuJDZ7fnVZrPMXW7cz0yAH5anYnWeF0fQNPrIlwxbmDGNhLf2fbQkHSwX3U4DlcU8vB4NRadDoaPI3T7x9uPnimDe/HlZMGc+WkwYwd1EtHMpKRXlGDxnaT1CAxs6uBfwVygH93979rsvx+4G6gDjgAfN7dd5vZdOCHQB+gHvhrd18SvGc08DhQAKwBPuPutcn8HJ1JW4NnU+UxXty0j4df2MrDL2xl9MB85gehMmNEf3K6KFQkM5SWhxnUuxtzxqtBY1uZuydnxWY5wDZgPhAByoFPuvvmmDFzgTfd/YSZfQm4zN0Xm9l4wN19u5kNAyqAc939AzMrBZ5w98fN7EfAenf/YUu1lJSUeCgUSsrnlPj2HT3JS1v289Lm/byx8yCn650B+Xlcce4grpw0hIvHDaR715x0lymd1P5jJ7nob1fyhTnn8PWrJ6a7nIxlZhXuXnK2cck8IpkN7HD3XUFBjwM3AmeCxN1XxYxfDdwezN8WM6bSzKqAQjM7ClwOfCpY/F/AXxI9epEMMqRvdz5z4Ug+c+FIjp08za+2HuDFzft57u19lIYi9Oiaw6XjBzJ/0hCumDiI/vl56S5ZOpFlayI0ODqt1U6SGSRFQDjmdQS4oIXxdwHPNZ1pZrOBPGAnMAD4wN0bT8xHgu1IBuvTvSufmDaMT0wbRm1dA6t3HeKlzdGjlRc27aeLwfmjCrhy8hCunDSY4QU9012ydGDRBo0RZo8uYPTA/HSX0yEkM0jinQyPex7NzG4HSoA5TeYPBX4O3OHuDRb/qm1z67wHuAdgxIgRrShbkikvtwuXji/k0vGFPHTjZN7ec5SXNu/nxU37+fYzm/n2M5uZOKQ3V04azPxJQzivqI8u1ku7Kn/vCO8erOErc8emu5QOI5lBEgFijxuLgcqmg8xsHvAgMMfdT8XM7wM8C3zD3VcHsw8C/cwsNzgqibtOAHd/BHgEotdI2v5xpL2ZGVOL+zG1uB//98oJ7D5UEw2Vzfv5t1U7+O4rOxjatzvzJw1m/qTBXDB6AHm5algtbbOkPEyvbrlcO2VIukvpMJIZJOXAuOAuqz3Abfzu2gYAZjYD+DFwtbtXxczPA54EfubuZY3z3d3NbBWwgOidW3cAv0ziZ5AUGjkgn7svGcPdl4zhcE0tK4OL9aWhMD97Yze9u+cyd8Igrpw8mDnjC+ndvWu6S5Ysc/zkaVa8vZebZqhBY3tK2p509zozuxd4gejtv4+6+yYzewgIufty4GGgF1AWnL54391vABYBlwIDzOzOYJV3uvs64OvA42b2/4G1wH8k6zNI+hTk57GwZDgLS4bzYW09v9lxkBc37WPlO1UsX19J1xzjonMGBqfABjO4T/d0lyxZ4JkNe6MNGnWRvV0l7fbfTKLbfzuO+ganYvcRXtq8jxc372f3oROAHoKUxNz8g9epPlnHi/ddqr8jCUj09l8FiWQtd2d7VXVwsX4f6yNHARg1oCdXTh7C/EmDmamHICWwff9x5v/za3zjunO5+5Ix6S4nK2TCcyQiSWVmjB/cm/GDe/OVuWN/7yHI/3z9XR55bdeZhyDnTxrCJXoIslMrDYXJ7WLcNENPDLQ3BYl0GPEegnypyUOQl4wbyJWTh3D5xEEU6CHITqO2roEn1uxh3rmD1aAxCRQk0iE1fQjyzXcP8eKm/WduL+5iUDKqILiuMoQRA/QQZEf2yjtVHKqpZdH5xekupUPSNRLpVNz9zEOQL23ezzv7jgMwcUjvoLmkHoLsiD7/03I2VR7l9a9fTm6OnkVKlK6RiMTR9CHI9w+d4MXgDrDvr9rB9/QQZIez7+hJXt1axRfnnKMQSRIFiXRqIwb0TOghyPmTBnPZBD0EmY3UoDH5FCQiAT0E2fFEGzSGuWB0AaPUoDFpFCQicfTIyzlzeqvpQ5DfeGoj33hqI9OK+555XmWcHoLMSG+9e5j3Dp3g/1w+Lt2ldGi62C7SCi09BDl/0mCunDxED0FmkPtL1/Hipv2UPziPHnl6hqi1dLFdJAlaegjyp799j5/8+l09BJkhGhs03jyjWCGSZAoSkTY420OQ3bt24dJxhcyfNJgrzh2shyBT6On1ezl5uoHF5+sie7IpSETaiR6CzCyloTDjB/diWnHfdJfS4ekaiUiSNfcQ5ITBvbnqvCHcOrOIkQN0R1F72rb/OFeqQWOb6RqJSIZo6SHI772yne+u3M7s0QUsmFXMdVOGkt9N/yzbqrQ8TNcc42Y1aEwJHZGIpNHeox/yxJo9LK2I8O7BGnrm5XDtlKEsmFXMBaMLdEvxR1Bb18CFf7uSC0YX8MPbZ6W7nKymIxKRLDC0bw++MncsX77sHCp2H6EsFOGZDZUsrYgwoqAnC2YVc8vMIor763pKolZu2c/hmlo9yZ5COiIRyTAnaut4fuM+ykIR3th1CDP42DkDWDhrOFdNHqJbWc/ic//5Flv2Huf1By7X8zxtpCMSkSzVMy+XW2YWc8vMYsKHT7BsTYSlFRG+tmQdvbvlcv20oSyYNZyZI/rp1FcT+46e5FfbDvDly8YqRFJIQSKSwYYX9ORr88bzx5eP4813D1NWEeaptZX8z1thxhTmR099zShmSF/1/YLfNWhcWKLvHUklndoSyTLVp+pYsWEvZRVhyt87QheDS8YVsrCkmHnnDu60T9I3NDhz//FVhvbtzuP3XJTucjoEndoS6aB6dctl0fnDWXT+cN47WMPSigjL1kS497/X0rdHV26YNoyFJcVMKerbqU59vfXeYXYfOsFXr1CDxlRTkIhksVED8/mTqyZw3/zx/HbnQcpCEUpDYX6+ejcTBvdmwaxibppRRGHvjv895aXlYXp3y+Wa84amu5ROR0Ei0gHkdDEuGVfIJeMKOfrhaZ7ZUElZKMJfr9jC3z3/DnMnFLJg1nAunzioQ37j47GTp1mxcS+3zFSDxnRQkIh0MH17dOXTF4zk0xeMZEfVccoqIjyxZg8vb6miID+PG6cPY+Gs4Uwa1ifdpbabp9dXRhs06tmRtNDFdpFOoK6+gV9vP0hZRZiXN1dRW9/A5GF9WDCrmBunF2V9V+Ibv/86J2vref5rl3Sq60LJpovtInJGbk4X5k4cxNyJgzhSU8vy9ZWUVYT5q6c38zcrtjDv3MEsmFXMnPGF5OZk16mvrfuOsz78AX9x/SSFSJooSEQ6mf75edzxsVHc8bFRbNl7jKUVEZ5au4fnNu6jsHc3bplRxIJZxYwb3DvdpSZkiRo0pp1ObYkItXUNrNpaxdKKCKveqaKuwZk2vB8LZhVzw9Rh9O3ZNd0lxlVb18AFf/MyF50zgB98Wg0a25tObYlIwvJyu3DV5CFcNXkIB6tP8dTaaEfiv3hqI99+ZjNXTR7CglnFXDx2YEa1Hnl5y36OnDjNQl1kTysFiYj8noG9unH3JWO46+LRbKo8RlkozC/XV/L0+kqG9u3OLTOLuHVmMWMKe6W7VEpDYYb27c6l4wrTXUqnpiARkbjMjPOK+nJeUV/+/LpzWbmlirJQmB++upPvr9pJycj+LCwp5topQ+ndPfWnvvYe/ZDXth3gK3PVoDHdFCQiclbdcqNfuHXtlKHsP3aSJ9fuoSwU5uvL3uYvl2/mmvOGsKCkmAtHD6BLin6oL6sIGjTO0mmtdFOQiEirDO7TnS/OOYcvXDqGteEPWFoR4el1lTyxdg/F/Xtw68xiFswqZnhB8r6Mq6HBKQ1FuGjMAEYM0Jd+pZuCREQ+EjNj5oj+zBzRn29eP4kXNu1jaUWE776ynX9duZ0LxxSwcNZwrpkyhJ557fuj5s13D/P+4RPcN18NGjOBbv8VkXa154MPeTL4Mq73Dp0gPy+H66YOZWHJcEpG9m+XhwbvW7KOl7fsp/zBeZ22bX4q6PZfEUmLon49uPfycXxl7ljK3zvC0oowz27YS2kowqgBjd9DX8ywfj0+0vqPnTzNirf3smBWsUIkQyhIRCQpzIzZowuYPbqAb31iMs9t3MfSijD/8OI2/vGlbVw8diALZhVz1eQhrQqE5esqOVXXwOLzdZE9UyhIRCTp8rvlsmBW9CL8+4dOsHRNhGUVEb76+Dp6d8/lhmnDWDCrmOnDz/499GWhMBOH9GZKUd8UVS9nk9TubGZ2tZltNbMdZvZAnOX3m9lmM9tgZivNbGTMsufN7AMze6bJe64wszVmts7MfmNmY5P5GUSkfY0Y0JP754/n1382l/+++wLmnTuYZWsi3PyD3zL/n1/jR7/aSdWxk3Hf+86+Y6yPHGVRyXA1aMwgSbvYbmY5wDZgPhAByoFPuvvmmDFzgTfd/YSZfQm4zN0XB8uuAHoCX3D362Pesw240d23mNmXgdnufmdLtehiu0hmO37yNM9u2MvSigih3UfI6WLMGV/IglnFXHHuILrlRk99/dXTm3hs9W7e/PN5Wd/6PhtkwsX22cAOd98VFPQ4cCNwJkjcfVXM+NXA7THLVprZZXHW60DjN/L0BSrbt2wRSbXe3bty2+wR3DZ7BLsOVLNsTYRlFXv48jtr6NezKzdNL+LG6cN4au0erpw0RCGSYZIZJEVAOOZ1BLighfF3Ac8lsN67gRVm9iFwDLgw3iAzuwe4B2DEiBGJ1CsiGWBMYS/+9KqJ3D9/Ar/ZcZClFRH++633+elv3wNgYUlxeguUP5DMIIl3AjPueTQzux3O9JYuAAAIoElEQVQoAeYksN77gGvd/U0z+1Pgn4iGy+9vyP0R4BGIntpKtGgRyQyNp7fmjC/k6InTLN9QyZ4jH3KJGjRmnGQGSQSIvT+vmDinocxsHvAgMMfdT7W0QjMrBKa5+5vBrCXA8+1Trohkqr49u/KZC0eefaCkRTLv2ioHxpnZaDPLA24DlscOMLMZwI+BG9y9KoF1HgH6mtn44PV8YEs71iwiIq2UtCMSd68zs3uBF4Ac4FF332RmDwEhd18OPAz0AsqCW/ned/cbAMzs18BEoJeZRYC73P0FM/sjYJmZNRANls8n6zOIiMjZqdeWiIjElejtv0l9IFFERDo+BYmIiLSJgkRERNpEQSIiIm2iIBERkTbpFHdtmdkBYPdHfPtA4GA7ltNeVFfrqK7WUV2t01HrGunuZ20l0CmCpC3MLJTI7W+pprpaR3W1jupqnc5el05tiYhImyhIRESkTRQkZ/dIugtohupqHdXVOqqrdTp1XbpGIiIibaIjEhERaRMFSRNm9rCZvWNmG8zsSTPr18y4q81sq5ntMLMHUlDXQjPbZGYNZtbsXRhm9p6ZvW1m68ws6Z0qW1FXqvdXgZm9ZGbbg//2b2ZcfbCv1pnZ8nhj2qmeFj+/mXUzsyXB8jfNbFSyamllXXea2YGYffQHXyKXpLoeNbMqM9vYzHIzs+8GdW8ws5kZUNNlZnY0Zl99M9k1BdsdbmarzGxL8G/xq3HGJHd/ubv+xPwBrgRyg+m/B/4+zpgcYCcwBsgD1gOTklzXucAE4FWgpIVx7wEDU7i/zlpXmvbXd4AHgukH4v1/DJZVp2AfnfXzA18GfhRM3wYsyZC67gT+LVV/n2K2eykwE9jYzPJriX41txH9uu03M6Cmy4Bn0rCvhgIzg+newLY4/x+Tur90RNKEu7/o7nXBy9VEv9mxqdnADnff5e61wOPAjUmua4u7b03mNj6KBOtK+f4K1v9fwfR/ATcleXstSeTzx9a7FLjCgi/pSXNdaeHurwGHWxhyI/Azj1oN9DOzoWmuKS3cfa+7rwmmjxP9sr+iJsOSur8UJC37PNEUb6oICMe8jvCH/+PSxYEXzazCzO5JdzGBdOyvwe6+F6L/0IBBzYzrbmYhM1ttZskKm0Q+/5kxwS8yR4EBSaqnNXUB3BqcDllqZsPjLE+HTP03eJGZrTez58xscqo3HpwSnQG82WRRUvdXMr+zPWOZ2cvAkDiLHnT3XwZjHgTqgF/EW0WceW2+/S2RuhLwcXevNLNBwEtm9k7wm1Q660r5/mrFakYE+2sM8IqZve3uO9taWxOJfP6k7KOzSGSbTwP/4+6nzOyLRI+aLk9yXYlIx/46mzVEW4pUm9m1wFPAuFRt3Mx6AcuAr7n7saaL47yl3fZXpwwSd5/X0nIzuwO4HrjCgxOMTUSA2N/MioHKZNeV4Doqg/9WmdmTRE9ftClI2qGulO8vM9tvZkPdfW9wCF/VzDoa99cuM3uV6G9z7R0kiXz+xjERM8sF+pL80yhnrcvdD8W8/AnR64aZICl/p9oi9oe3u68wsx+Y2UB3T3oPLjPrSjREfuHuT8QZktT9pVNbTZjZ1cDXgRvc/UQzw8qBcWY22szyiF4cTdodP4kys3wz6904TfTGgbh3mKRYOvbXcuCOYPoO4A+OnMysv5l1C6YHAh8HNiehlkQ+f2y9C4BXmvklJqV1NTmPfgPR8++ZYDnw2eBupAuBo42nMtPFzIY0Xtcys9lEf74eavld7bJdA/4D2OLu/9TMsOTur1TfYZDpf4AdRM8lrgv+NN5JMwxYETPuWqJ3R+wkeoon2XXdTPS3ilPAfuCFpnURvftmffBnU6bUlab9NQBYCWwP/lsQzC8B/j2Y/hjwdrC/3gbuSmI9f/D5gYeI/sIC0B0oC/7+vQWMSfY+SrCuvw3+Lq0HVgETU1TX/wB7gdPB36+7gC8CXwyWG/D9oO63aeFOxhTWdG/MvloNfCxF++pioqepNsT83Lo2lftLT7aLiEib6NSWiIi0iYJERETaREEiIiJtoiAREZE2UZCIiEibKEhEWmBm1W18/9LgqfmWxrxqLXROTnRMk/GFZvZ8ouNF2kJBIpIkQa+lHHffleptu/sBYK+ZfTzV25bOR0EikoDgieCHzWyjRb/vZXEwv0vQCmOTmT1jZivMbEHwtk8T80S9mf0waBC5ycz+qpntVJvZP5rZGjNbaWaFMYsXmtlbZrbNzC4Jxo8ys18H49eY2cdixj8V1CCSVAoSkcTcAkwHpgHzgIeD9iG3AKOAKcDdwEUx7/k4UBHz+kF3LwGmAnPMbGqc7eQDa9x9JvAr4Fsxy3LdfTbwtZj5VcD8YPxi4Lsx40PAJa3/qCKt0ymbNop8BBcT7YJbD+w3s18B5wfzy9y9AdhnZqti3jMUOBDzelHQ2j83WDaJaFuLWA3AkmD6MSC2AV/jdAXR8ALoCvybmU0H6oHxMeOriLaqEUkqBYlIYpr7kqmWvnzqQ6I9tDCz0cCfAOe7+xEz+2njsrOI7WF0KvhvPb/7t3sf0R5n04ieYTgZM757UINIUunUlkhiXgMWm1lOcN3iUqLNFX9D9IufupjZYKJft9poCzA2mO4D1ABHg3HXNLOdLkS7/wJ8Klh/S/oCe4Mjos8Q/frcRuPJjO7P0sHpiEQkMU8Svf6xnuhRwp+5+z4zWwZcQfQH9jai30x3NHjPs0SD5WV3X29ma4l2h90FvN7MdmqAyWZWEaxn8Vnq+gGwzMwWEu3OWxOzbG5Qg0hSqfuvSBuZWS+PfiveAKJHKR8PQqYH0R/uHw+urSSyrmp379VOdb0G3OjuR9pjfSLN0RGJSNs9Y2b9gDzg2+6+D8DdPzSzbxH9buz3U1lQcPrtnxQikgo6IhERkTbRxXYREWkTBYmIiLSJgkRERNpEQSIiIm2iIBERkTZRkIiISJv8L5TEVN3D8AwSAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "('alpha is:', 10.0)\n"
     ]
    }
   ],
   "source": [
    "mse_mean = np.mean(ridge.cv_values_, axis = 0)\n",
    "plt.plot(np.log10(alphas), mse_mean.reshape(len(alphas),1)) \n",
    "\n",
    "#这是为了标出最佳参数的位置，不是必须\n",
    "#plt.plot(np.log10(ridge.alpha_)*np.ones(3), [0.28, 0.29, 0.30])\n",
    "\n",
    "plt.xlabel('log(alpha)')\n",
    "plt.ylabel('mse')\n",
    "plt.show()\n",
    "\n",
    "print ('alpha is:', ridge.alpha_)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/anaconda2/lib/python2.7/site-packages/sklearn/linear_model/coordinate_descent.py:1094: DataConversionWarning: A column-vector y was passed when a 1d array was expected. Please change the shape of y to (n_samples, ), for example using ravel().\n",
      "  y = column_or_1d(y, warn=True)\n",
      "/anaconda2/lib/python2.7/site-packages/sklearn/linear_model/coordinate_descent.py:491: ConvergenceWarning: Objective did not converge. You might want to increase the number of iterations. Fitting data with very small alpha may cause precision problems.\n",
      "  ConvergenceWarning)\n",
      "/anaconda2/lib/python2.7/site-packages/sklearn/linear_model/coordinate_descent.py:491: ConvergenceWarning: Objective did not converge. You might want to increase the number of iterations. Fitting data with very small alpha may cause precision problems.\n",
      "  ConvergenceWarning)\n",
      "/anaconda2/lib/python2.7/site-packages/sklearn/linear_model/coordinate_descent.py:491: ConvergenceWarning: Objective did not converge. You might want to increase the number of iterations. Fitting data with very small alpha may cause precision problems.\n",
      "  ConvergenceWarning)\n",
      "/anaconda2/lib/python2.7/site-packages/sklearn/linear_model/coordinate_descent.py:491: ConvergenceWarning: Objective did not converge. You might want to increase the number of iterations. Fitting data with very small alpha may cause precision problems.\n",
      "  ConvergenceWarning)\n"
     ]
    }
   ],
   "source": [
    "#### Lasso／L1正则\n",
    "#class sklearn.linear_model.LassoCV(eps=0.001, n_alphas=100, alphas=None, fit_intercept=True, \n",
    " #                                   normalize=False, precompute=’auto’, max_iter=1000, \n",
    "  #                                  tol=0.0001, copy_X=True, cv=None, verbose=False, n_jobs=1,\n",
    "   #                                 positive=False, random_state=None, selection=’cyclic’)\n",
    "from sklearn.linear_model import LassoCV\n",
    "\n",
    "#设置超参数搜索范围\n",
    "#alphas = [ 0.01, 0.1, 1, 10,100]\n",
    "\n",
    "#生成一个LassoCV实例\n",
    "#lasso = LassoCV(alphas=alphas)  \n",
    "lasso = LassoCV()  \n",
    "\n",
    "#训练（内含CV）\n",
    "lasso.fit(X_train, y_train)  \n",
    "\n",
    "#测试\n",
    "y_test_pred_lasso = lasso.predict(X_test)\n",
    "y_train_pred_lasso = lasso.predict(X_train)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "('MSE:', 0.18149727035295665)\n",
      "('RMSE:', 0.4260249644715162)\n"
     ]
    }
   ],
   "source": [
    "print (\"MSE:\",metrics.mean_squared_error(y_test, y_test_pred_lasso))\n",
    "\n",
    "print (\"RMSE:\",np.sqrt(metrics.mean_squared_error(y_test, y_test_pred_lasso)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "从以上实验的三个方法中可以看出岭回归模型的RMSE的值最小，说明岭回归模型对于这个例子效果最好。岭回归的目标函数是在OLS上加了一个正则项."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.15"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
